diff --git a/NPLib/Core/NPOptionManager.cxx b/NPLib/Core/NPOptionManager.cxx
index 79f648d3175074a53a7ef494a4d04585c7f03c41..c4e1573dcac0485d213d57210bd726d570ced8c7 100644
--- a/NPLib/Core/NPOptionManager.cxx
+++ b/NPLib/Core/NPOptionManager.cxx
@@ -1,4 +1,4 @@
-/*****************************************************************************
+;/*****************************************************************************
* Copyright (C) 2009-2016 this file is part of the NPTool Project *
* *
* For the licensing terms see $NPTOOL/Licence/NPTool_Licence *
@@ -96,7 +96,7 @@ void NPOptionManager::ReadTheInputArgument(int argc, char** argv){
- for (int i = 0; i < argc; i++) {
+ for (int i = 1; i < argc; i++) {
std::string argument = argv[i];
if (argument == "-H" || argument == "-h" || argument == "--help") DisplayHelp();
@@ -168,8 +168,9 @@ void NPOptionManager::ReadTheInputArgument(int argc, char** argv){
else if (argument == "--circular") {fCircularTree = true;}
-
- //else ;
+ else{
+ SendErrorAndExit(argument.c_str());
+ }
}
CheckArguments();
if(argc!=0)
@@ -384,7 +385,8 @@ void NPOptionManager::SendErrorAndExit(const char* type) const{
}
else {
- std::cout << "NPOptionManager::SendErrorAndExit() unkwown keyword" << std::endl;
+ std::cout << "NPOptionManager::SendErrorAndExit() unkwown program argument: " << type << std::endl;
+ exit(1);
}
}
diff --git a/NPLib/Detectors/PISTA/CMakeLists.txt b/NPLib/Detectors/PISTA/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..fe7686eee6bfe0d91d709d23238b4c5153050672
--- /dev/null
+++ b/NPLib/Detectors/PISTA/CMakeLists.txt
@@ -0,0 +1,6 @@
+add_custom_command(OUTPUT TPISTAPhysicsDict.cxx COMMAND ../../scripts/build_dict.sh TPISTAPhysics.h TPISTAPhysicsDict.cxx TPISTAPhysics.rootmap libNPPISTA.dylib DEPENDS TPISTAPhysics.h)
+add_custom_command(OUTPUT TPISTADataDict.cxx COMMAND ../../scripts/build_dict.sh TPISTAData.h TPISTADataDict.cxx TPISTAData.rootmap libNPPISTA.dylib DEPENDS TPISTAData.h)
+add_library(NPPISTA SHARED TPISTASpectra.cxx TPISTAData.cxx TPISTAPhysics.cxx TPISTADataDict.cxx TPISTAPhysicsDict.cxx )
+target_link_libraries(NPPISTA ${ROOT_LIBRARIES} NPCore)
+install(FILES TPISTAData.h TPISTAPhysics.h TPISTASpectra.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY})
+
diff --git a/NPLib/Detectors/PISTA/TPISTAData.cxx b/NPLib/Detectors/PISTA/TPISTAData.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..5a5970934c58048f7570b77e68b6693fcf7898f4
--- /dev/null
+++ b/NPLib/Detectors/PISTA/TPISTAData.cxx
@@ -0,0 +1,107 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold PISTA Raw data *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+#include "TPISTAData.h"
+
+#include <iostream>
+#include <fstream>
+#include <sstream>
+#include <string>
+using namespace std;
+
+ClassImp(TPISTAData)
+
+
+//////////////////////////////////////////////////////////////////////
+TPISTAData::TPISTAData() {
+}
+
+
+
+//////////////////////////////////////////////////////////////////////
+TPISTAData::~TPISTAData() {
+}
+
+
+
+//////////////////////////////////////////////////////////////////////
+void TPISTAData::Clear() {
+ // Energy X
+ fFirstStage_XE_DetectorNbr.clear();
+ fFirstStage_XE_StripNbr.clear();
+ fFirstStage_XE_Energy.clear();
+ // Energy Y
+ fFirstStage_YE_DetectorNbr.clear();
+ fFirstStage_YE_StripNbr.clear();
+ fFirstStage_YE_Energy.clear();
+ // Time X
+ fFirstStage_XT_DetectorNbr.clear();
+ fFirstStage_XT_StripNbr.clear();
+ fFirstStage_XT_Time.clear();
+ // Time Y
+ fFirstStage_YT_DetectorNbr.clear();
+ fFirstStage_YT_StripNbr.clear();
+ fFirstStage_YT_Time.clear();
+
+ // Energy X
+ fSecondStage_XE_DetectorNbr.clear();
+ fSecondStage_XE_StripNbr.clear();
+ fSecondStage_XE_Energy.clear();
+ // Energy Y
+ fSecondStage_YE_DetectorNbr.clear();
+ fSecondStage_YE_StripNbr.clear();
+ fSecondStage_YE_Energy.clear();
+ // Time X
+ fSecondStage_XT_DetectorNbr.clear();
+ fSecondStage_XT_StripNbr.clear();
+ fSecondStage_XT_Time.clear();
+ // Time Y
+ fSecondStage_YT_DetectorNbr.clear();
+ fSecondStage_YT_StripNbr.clear();
+ fSecondStage_YT_Time.clear();
+
+}
+
+
+
+//////////////////////////////////////////////////////////////////////
+void TPISTAData::Dump() const {
+ // This method is very useful for debuging and worth the dev.
+ cout << "XXXXXXXXXXXXXXXXXXXXXXXX New Event [TPISTAData::Dump()] XXXXXXXXXXXXXXXXX" << endl;
+
+ // Energy
+ size_t mysize = fFirstStage_XE_DetectorNbr.size();
+ cout << "First Stage PISTA_XE_Mult: " << mysize << endl;
+
+ for (size_t i = 0 ; i < mysize ; i++){
+ cout << "X-DetNbr: " << fFirstStage_XE_DetectorNbr[i]
+ << " X-Energy: " << fFirstStage_XE_Energy[i];
+ }
+
+ // Time
+ mysize = fFirstStage_XT_DetectorNbr.size();
+ cout << "First Stage PISTA_XT_Mult: " << mysize << endl;
+
+ for (size_t i = 0 ; i < mysize ; i++){
+ cout << "X-DetNbr: " << fFirstStage_XT_DetectorNbr[i]
+ << " X-Time: " << fFirstStage_XT_Time[i];
+ }
+}
diff --git a/NPLib/Detectors/PISTA/TPISTAData.h b/NPLib/Detectors/PISTA/TPISTAData.h
new file mode 100644
index 0000000000000000000000000000000000000000..82f0a3861ec51b516bea5432d803167448709eaf
--- /dev/null
+++ b/NPLib/Detectors/PISTA/TPISTAData.h
@@ -0,0 +1,227 @@
+#ifndef __PISTADATA__
+#define __PISTADATA__
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold PISTA Raw data *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+// STL
+#include <vector>
+using namespace std;
+
+// ROOT
+#include "TObject.h"
+
+class TPISTAData : public TObject {
+ //////////////////////////////////////////////////////////////
+ // data members are hold into vectors in order
+ // to allow multiplicity treatment
+ private:
+ // First Stage Front Energy
+ vector<unsigned short> fFirstStage_XE_DetectorNbr;
+ vector<unsigned short> fFirstStage_XE_StripNbr;
+ vector<double> fFirstStage_XE_Energy;
+ // First Stage Front Time
+ vector<unsigned short> fFirstStage_XT_DetectorNbr;
+ vector<unsigned short> fFirstStage_XT_StripNbr;
+ vector<double> fFirstStage_XT_Time;
+ // First Stage Back Energy
+ vector<unsigned short> fFirstStage_YE_DetectorNbr;
+ vector<unsigned short> fFirstStage_YE_StripNbr;
+ vector<double> fFirstStage_YE_Energy;
+ // First Stage Back Time
+ vector<unsigned short> fFirstStage_YT_DetectorNbr;
+ vector<unsigned short> fFirstStage_YT_StripNbr;
+ vector<double> fFirstStage_YT_Time;
+
+ // Second Stage Front Energy
+ vector<unsigned short> fSecondStage_XE_DetectorNbr;
+ vector<unsigned short> fSecondStage_XE_StripNbr;
+ vector<double> fSecondStage_XE_Energy;
+ // Second Stage Front Time
+ vector<unsigned short> fSecondStage_XT_DetectorNbr;
+ vector<unsigned short> fSecondStage_XT_StripNbr;
+ vector<double> fSecondStage_XT_Time;
+ // Second Stage Back Energy
+ vector<unsigned short> fSecondStage_YE_DetectorNbr;
+ vector<unsigned short> fSecondStage_YE_StripNbr;
+ vector<double> fSecondStage_YE_Energy;
+ // Second Stage Back Time
+ vector<unsigned short> fSecondStage_YT_DetectorNbr;
+ vector<unsigned short> fSecondStage_YT_StripNbr;
+ vector<double> fSecondStage_YT_Time;
+
+
+ //////////////////////////////////////////////////////////////
+ // Constructor and destructor
+ public:
+ TPISTAData();
+ ~TPISTAData();
+
+
+ //////////////////////////////////////////////////////////////
+ // Inherited from TObject and overriden to avoid warnings
+ public:
+ void Clear();
+ void Clear(const Option_t*) {};
+ void Dump() const;
+
+
+ //////////////////////////////////////////////////////////////
+ // Getters and Setters
+ // Prefer inline declaration to avoid unnecessary called of
+ // frequently used methods
+ // add //! to avoid ROOT creating dictionnary for the methods
+ public:
+ ////////////////////// SETTERS ////////////////////////
+ // First Stage Energy Front
+ inline void SetFirstStageXE(const UShort_t& DetNbr, const UShort_t& StripNbr, const Double_t& Energy){
+ fFirstStage_XE_DetectorNbr.push_back(DetNbr);
+ fFirstStage_XE_StripNbr.push_back(StripNbr);
+ fFirstStage_XE_Energy.push_back(Energy);
+ };//!
+ // First Stage Energy Back
+ inline void SetFirstStageYE(const UShort_t& DetNbr, const UShort_t& StripNbr, const Double_t& Energy){
+ fFirstStage_YE_DetectorNbr.push_back(DetNbr);
+ fFirstStage_YE_StripNbr.push_back(StripNbr);
+ fFirstStage_YE_Energy.push_back(Energy);
+ };//!
+ // First Stage Time Front
+ inline void SetFirstStageXT(const UShort_t& DetNbr, const UShort_t StripNbr, const Double_t& Time) {
+ fFirstStage_XT_DetectorNbr.push_back(DetNbr);
+ fFirstStage_XT_StripNbr.push_back(StripNbr);
+ fFirstStage_XT_Time.push_back(Time);
+ };//!
+ // First Stage Time Back
+ inline void SetFirstStageYT(const UShort_t& DetNbr, const UShort_t StripNbr, const Double_t& Time) {
+ fFirstStage_YT_DetectorNbr.push_back(DetNbr);
+ fFirstStage_YT_StripNbr.push_back(StripNbr);
+ fFirstStage_YT_Time.push_back(Time);
+ };//!
+
+ //////
+ // Second Stage Energy Front
+ inline void SetSecondStageXE(const UShort_t& DetNbr, const UShort_t& StripNbr, const Double_t& Energy){
+ fSecondStage_XE_DetectorNbr.push_back(DetNbr);
+ fSecondStage_XE_StripNbr.push_back(StripNbr);
+ fSecondStage_XE_Energy.push_back(Energy);
+ };//!
+ // Second Stage Energy Back
+ inline void SetSecondStageYE(const UShort_t& DetNbr, const UShort_t& StripNbr, const Double_t& Energy){
+ fSecondStage_YE_DetectorNbr.push_back(DetNbr);
+ fSecondStage_YE_StripNbr.push_back(StripNbr);
+ fSecondStage_YE_Energy.push_back(Energy);
+ };//!
+ // Second Stage Time Front
+ inline void SetSecondStageXT(const UShort_t& DetNbr, const UShort_t StripNbr, const Double_t& Time) {
+ fSecondStage_XT_DetectorNbr.push_back(DetNbr);
+ fSecondStage_XT_StripNbr.push_back(StripNbr);
+ fSecondStage_XT_Time.push_back(Time);
+ };//!
+ // Second Stage Time Back
+ inline void SetSecondStageYT(const UShort_t& DetNbr, const UShort_t StripNbr, const Double_t& Time) {
+ fSecondStage_YT_DetectorNbr.push_back(DetNbr);
+ fSecondStage_YT_StripNbr.push_back(StripNbr);
+ fSecondStage_YT_Time.push_back(Time);
+ };//!
+
+ ////////////////////// GETTERS ////////////////////////
+ // First Stage Energy X
+ inline UShort_t GetFirstStageMultXEnergy() const
+ {return fFirstStage_XE_DetectorNbr.size();}
+ inline UShort_t GetFirstStage_XE_DetectorNbr(const unsigned int &i) const
+ {return fFirstStage_XE_DetectorNbr[i];}//!
+ inline UShort_t GetFirstStage_XE_StripNbr(const unsigned int &i) const
+ {return fFirstStage_XE_StripNbr[i];}//!
+ inline Double_t GetFirstStage_XE_Energy(const unsigned int &i) const
+ {return fFirstStage_XE_Energy[i];}//!
+ // First Stage Energy Y
+ inline UShort_t GetFirstStageMultYEnergy() const
+ {return fFirstStage_YE_DetectorNbr.size();}
+ inline UShort_t GetFirstStage_YE_DetectorNbr(const unsigned int &i) const
+ {return fFirstStage_YE_DetectorNbr[i];}//!
+ inline UShort_t GetFirstStage_YE_StripNbr(const unsigned int &i) const
+ {return fFirstStage_YE_StripNbr[i];}//!
+ inline Double_t GetFirstStage_YE_Energy(const unsigned int &i) const
+ {return fFirstStage_YE_Energy[i];}//!
+ // First Stage Time X
+ inline UShort_t GetFirstStageMultXTime() const
+ {return fFirstStage_XT_DetectorNbr.size();}
+ inline UShort_t GetFirstStage_XT_DetectorNbr(const unsigned int &i) const
+ {return fFirstStage_XT_DetectorNbr[i];}//!
+ inline UShort_t GetFirstStage_XT_StripNbr(const unsigned int &i) const
+ {return fFirstStage_XT_StripNbr[i];}//!
+ inline Double_t GetFirstStage_XT_Time(const unsigned int &i) const
+ {return fFirstStage_XT_Time[i];}//!
+ // First Stage Time Y
+ inline UShort_t GetFirstStageMultYTime() const
+ {return fFirstStage_YT_DetectorNbr.size();}
+ inline UShort_t GetFirstStage_YT_DetectorNbr(const unsigned int &i) const
+ {return fFirstStage_YT_DetectorNbr[i];}//!
+ inline UShort_t GetFirstStage_YT_StripNbr(const unsigned int &i) const
+ {return fFirstStage_YT_StripNbr[i];}//!
+ inline Double_t GetFirstStage_YT_Time(const unsigned int &i) const
+ {return fFirstStage_YT_Time[i];}//!
+
+ //////
+ // Second Stage Energy X
+ inline UShort_t GetSecondStageMultXEnergy() const
+ {return fSecondStage_XE_DetectorNbr.size();}
+ inline UShort_t GetSecondStage_XE_DetectorNbr(const unsigned int &i) const
+ {return fSecondStage_XE_DetectorNbr[i];}//!
+ inline UShort_t GetSecondStage_XE_StripNbr(const unsigned int &i) const
+ {return fSecondStage_XE_StripNbr[i];}//!
+ inline Double_t GetSecondStage_XE_Energy(const unsigned int &i) const
+ {return fSecondStage_XE_Energy[i];}//!
+ // Second Stage Energy Y
+ inline UShort_t GetSecondStageMultYEnergy() const
+ {return fSecondStage_YE_DetectorNbr.size();}
+ inline UShort_t GetSecondStage_YE_DetectorNbr(const unsigned int &i) const
+ {return fSecondStage_YE_DetectorNbr[i];}//!
+ inline UShort_t GetSecondStage_YE_StripNbr(const unsigned int &i) const
+ {return fSecondStage_YE_StripNbr[i];}//!
+ inline Double_t GetSecondStage_YE_Energy(const unsigned int &i) const
+ {return fSecondStage_YE_Energy[i];}//!
+ // Second Stage Time X
+ inline UShort_t GetSecondStageMultXTime() const
+ {return fSecondStage_XT_DetectorNbr.size();}
+ inline UShort_t GetSecondStage_XT_DetectorNbr(const unsigned int &i) const
+ {return fSecondStage_XT_DetectorNbr[i];}//!
+ inline UShort_t GetSecondStage_XT_StripNbr(const unsigned int &i) const
+ {return fSecondStage_XT_StripNbr[i];}//!
+ inline Double_t GetSecondStage_XT_Time(const unsigned int &i) const
+ {return fSecondStage_XT_Time[i];}//!
+ // Second Stage Time Y
+ inline UShort_t GetSecondStageMultYTime() const
+ {return fSecondStage_YT_DetectorNbr.size();}
+ inline UShort_t GetSecondStage_YT_DetectorNbr(const unsigned int &i) const
+ {return fSecondStage_YT_DetectorNbr[i];}//!
+ inline UShort_t GetSecondStage_YT_StripNbr(const unsigned int &i) const
+ {return fSecondStage_YT_StripNbr[i];}//!
+ inline Double_t GetSecondStage_YT_Time(const unsigned int &i) const
+ {return fSecondStage_YT_Time[i];}//!
+
+
+ //////////////////////////////////////////////////////////////
+ // Required for ROOT dictionnary
+ ClassDef(TPISTAData,1) // PISTAData structure
+};
+
+#endif
diff --git a/NPLib/Detectors/PISTA/TPISTAPhysics.cxx b/NPLib/Detectors/PISTA/TPISTAPhysics.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..7266d3d536f5b3d8325076a69385e997e9362e56
--- /dev/null
+++ b/NPLib/Detectors/PISTA/TPISTAPhysics.cxx
@@ -0,0 +1,552 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold PISTA Treated data *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+#include "TPISTAPhysics.h"
+
+// STL
+#include <sstream>
+#include <iostream>
+#include <cmath>
+#include <stdlib.h>
+#include <limits>
+using namespace std;
+
+// NPL
+#include "RootInput.h"
+#include "RootOutput.h"
+#include "NPDetectorFactory.h"
+#include "NPOptionManager.h"
+
+// ROOT
+#include "TChain.h"
+
+ClassImp(TPISTAPhysics)
+ ///////////////////////////////////////////////////////////////////////////
+ TPISTAPhysics::TPISTAPhysics(){
+ EventMultiplicity = 0;
+ m_EventData = new TPISTAData;
+ m_PreTreatedData = new TPISTAData;
+ m_EventPhysics = this;
+ m_Spectra = NULL;
+ m_E_RAW_Threshold = 0; // adc channels
+ m_E_Threshold = 0; // MeV
+ m_NumberOfDetectors = 0;
+ m_MaximumStripMultiplicityAllowed = 10;
+ m_StripEnergyMatching = 0.050;
+ }
+
+///////////////////////////////////////////////////////////////////////////
+/// A usefull method to bundle all operation to add a detector
+void TPISTAPhysics::AddDetector(TVector3){
+ // In That simple case nothing is done
+ // Typically for more complex detector one would calculate the relevant
+ // positions (stripped silicon) or angles (gamma array)
+ m_NumberOfDetectors++;
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::AddDetector(double R, double Theta, double Phi){
+ m_NumberOfDetectors++;
+
+ // Vector U on detector face (parallel to Y strips) Y strips are along X axis
+ TVector3 U;
+ // Vector V on detector face (parallel to X strips)
+ TVector3 V;
+ // Vector W normal to detector face (pointing to the back)
+ TVector3 W;
+ // Vector C position of detector face center
+ TVector3 C;
+
+ C = TVector3(R*sin(Theta)*cos(Phi),
+ R*sin(Theta)*sin(Phi),
+ R*cos(Theta));
+
+ TVector3 P = TVector3(cos(Theta)*cos(Phi),
+ cos(Theta)*sin(Phi),
+ -sin(Theta));
+
+ W = C.Unit();
+ U = W.Cross(P);
+ V = W.Cross(U);
+
+ U = U.Unit();
+ V = V.Unit();
+
+ double Height = 118; // mm
+ double Base = 95; // mm
+ double NumberOfStrips = 128;
+ double StripPitchHeight = Height / NumberOfStrips; // mm
+ double StripPitchBase = Base / NumberOfStrips; // mm
+
+ vector<double> lineX;
+ vector<double> lineY;
+ vector<double> lineZ;
+
+ vector<vector<double>> OneDetectorStripPositionX;
+ vector<vector<double>> OneDetectorStripPositionY;
+ vector<vector<double>> OneDetectorStripPositionZ;
+
+ double X, Y, Z;
+
+ // Moving C to the 1.1 Corner;
+ TVector3 Strip_1_1;
+ Strip_1_1 = C - (0.5*Base*U + 0.5*Height*V) + U*(StripPitchBase / 2.) + V*(StripPitchHeight / 2.);
+
+ TVector3 StripPos;
+ for(int i=0; i<NumberOfStrips; i++){
+ lineX.clear();
+ lineY.clear();
+ lineZ.clear();
+ for(int j=0; j<NumberOfStrips; j++){
+ StripPos = Strip_1_1 + i*U*StripPitchBase + j*V*StripPitchHeight;
+ lineX.push_back(StripPos.X());
+ lineY.push_back(StripPos.Y());
+ lineZ.push_back(StripPos.Z());
+ }
+
+ OneDetectorStripPositionX.push_back(lineX);
+ OneDetectorStripPositionY.push_back(lineY);
+ OneDetectorStripPositionZ.push_back(lineZ);
+ }
+
+ m_StripPositionX.push_back(OneDetectorStripPositionX);
+ m_StripPositionY.push_back(OneDetectorStripPositionY);
+ m_StripPositionZ.push_back(OneDetectorStripPositionZ);
+}
+
+///////////////////////////////////////////////////////////////////////////
+TVector3 TPISTAPhysics::GetPositionOfInteraction(const int i){
+ TVector3 Position = TVector3(GetStripPositionX(DetectorNumber[i], StripX[i], StripY[i]),
+ GetStripPositionY(DetectorNumber[i], StripX[i], StripY[i]),
+ GetStripPositionZ(DetectorNumber[i], StripX[i], StripY[i]));
+
+ return Position;
+}
+
+///////////////////////////////////////////////////////////////////////////
+TVector3 TPISTAPhysics::GetDetectorNormal(const int i){
+ TVector3 U = TVector3(GetStripPositionX(DetectorNumber[i],128,1),
+ GetStripPositionY(DetectorNumber[i],128,1),
+ GetStripPositionZ(DetectorNumber[i],128,1))
+
+ -TVector3(GetStripPositionX(DetectorNumber[i],1,1),
+ GetStripPositionY(DetectorNumber[i],1,1),
+ GetStripPositionZ(DetectorNumber[i],1,1));
+
+ TVector3 V = TVector3(GetStripPositionX(DetectorNumber[i],128,128),
+ GetStripPositionY(DetectorNumber[i],128,128),
+ GetStripPositionZ(DetectorNumber[i],128,128))
+
+ -TVector3(GetStripPositionX(DetectorNumber[i],128,1),
+ GetStripPositionY(DetectorNumber[i],128,1),
+ GetStripPositionZ(DetectorNumber[i],128,1));
+
+ TVector3 Normal = U.Cross(V);
+
+ return (Normal.Unit());
+}
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::BuildSimplePhysicalEvent() {
+ BuildPhysicalEvent();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::BuildPhysicalEvent() {
+ // apply thresholds and calibration
+ PreTreat();
+
+ if(1 /*CheckEvent() == 1*/){
+ vector<TVector2> couple = Match_X_Y();
+
+ EventMultiplicity = couple.size();
+ for(unsigned int i=0; i<couple.size(); i++){
+ int N = m_PreTreatedData->GetFirstStage_XE_DetectorNbr(couple[i].X());
+ int X = m_PreTreatedData->GetFirstStage_XE_StripNbr(couple[i].X());
+ int Y = m_PreTreatedData->GetFirstStage_YE_StripNbr(couple[i].Y());
+
+ double XE = m_PreTreatedData->GetFirstStage_XE_Energy(couple[i].X());
+ double YE = m_PreTreatedData->GetFirstStage_YE_Energy(couple[i].Y());
+ DetectorNumber.push_back(N);
+ StripX.push_back(X);
+ StripY.push_back(Y);
+ DE.push_back(XE);
+
+ PosX.push_back(GetPositionOfInteraction(i).x());
+ PosY.push_back(GetPositionOfInteraction(i).y());
+ PosZ.push_back(GetPositionOfInteraction(i).z());
+
+ int SecondStageMult = m_PreTreatedData->GetSecondStageMultXEnergy();
+ for(unsigned int j=0; j<SecondStageMult; j++){
+ if(m_PreTreatedData->GetSecondStage_XE_DetectorNbr(j)==N){
+ double XDE = m_PreTreatedData->GetSecondStage_XE_Energy(j);
+ double YDE = m_PreTreatedData->GetSecondStage_YE_Energy(j);
+
+ E.push_back(XDE);
+ }
+ }
+ }
+ }
+}
+
+///////////////////////////////////////////////////////////////////////////
+vector<TVector2> TPISTAPhysics::Match_X_Y(){
+ vector<TVector2> ArrayOfGoodCouple;
+
+ static unsigned int m_XEMult, m_YEMult;
+ m_XEMult = m_PreTreatedData->GetFirstStageMultXEnergy();
+ m_YEMult = m_PreTreatedData->GetFirstStageMultYEnergy();
+
+ if(m_XEMult>m_MaximumStripMultiplicityAllowed || m_YEMult>m_MaximumStripMultiplicityAllowed){
+ return ArrayOfGoodCouple;
+ }
+
+ for(unsigned int i=0; i<m_XEMult; i++){
+ for(unsigned int j=0; j<m_YEMult; j++){
+
+ // Declaration of variable for clarity
+ int XDetNbr = m_PreTreatedData->GetFirstStage_XE_DetectorNbr(i);
+ int YDetNbr = m_PreTreatedData->GetFirstStage_YE_DetectorNbr(j);
+
+ // if same detector check energy
+ if(XDetNbr == YDetNbr){
+ // Declaration of variable for clarity
+ double XE = m_PreTreatedData->GetFirstStage_XE_Energy(i);
+ double YE = m_PreTreatedData->GetFirstStage_YE_Energy(i);
+ double XStripNbr = m_PreTreatedData->GetFirstStage_XE_StripNbr(i);
+ double YStripNbr = m_PreTreatedData->GetFirstStage_YE_StripNbr(i);
+
+ // look if energy matches
+ if(abs(XE-YE)/2.<m_StripEnergyMatching){
+ ArrayOfGoodCouple.push_back(TVector2(i,j));
+ }
+ }
+ }
+ }
+
+ return ArrayOfGoodCouple;
+}
+
+///////////////////////////////////////////////////////////////////////////
+int TPISTAPhysics::CheckEvent(){
+ // Check the size of the different elements
+ if(m_PreTreatedData->GetFirstStageMultXEnergy() == m_PreTreatedData->GetFirstStageMultYEnergy() )
+ return 1;
+
+ else
+ return -1;
+}
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::PreTreat() {
+ // This method typically applies thresholds and calibrations
+ // Might test for disabled channels for more complex detector
+
+ // clear pre-treated object
+ ClearPreTreatedData();
+
+ // instantiate CalibrationManager
+ static CalibrationManager* Cal = CalibrationManager::getInstance();
+
+ //////
+ // First Stage Energy
+ unsigned int sizeFront = m_EventData->GetFirstStageMultXEnergy();
+ for (UShort_t i = 0; i < sizeFront ; ++i) {
+ if (m_EventData->GetFirstStage_XE_Energy(i) > m_E_RAW_Threshold) {
+ Double_t Energy = m_EventData->GetFirstStage_XE_Energy(i);
+ //Double_t Energy = Cal->ApplyCalibration("PISTA/ENERGY"+NPL::itoa(m_EventData->GetFirstStage_XE_DetectorNbr(i)),m_EventData->GetFirstStage_XE_Energy(i));
+ if (Energy > m_E_Threshold) {
+ m_PreTreatedData->SetFirstStageXE(m_EventData->GetFirstStage_XE_DetectorNbr(i), m_EventData->GetFirstStage_XE_StripNbr(i), Energy);
+ }
+ }
+ }
+ unsigned int sizeBack = m_EventData->GetFirstStageMultXEnergy();
+ for (UShort_t i = 0; i < sizeBack ; ++i) {
+ if (m_EventData->GetFirstStage_YE_Energy(i) > m_E_RAW_Threshold) {
+ Double_t Energy = m_EventData->GetFirstStage_YE_Energy(i);
+ //Double_t Energy = Cal->ApplyCalibration("PISTA/ENERGY"+NPL::itoa(m_EventData->GetFirstStage_YE_DetectorNbr(i)),m_EventData->GetFirstStage_YE_Energy(i));
+ if (Energy > m_E_Threshold) {
+ m_PreTreatedData->SetFirstStageYE(m_EventData->GetFirstStage_YE_DetectorNbr(i), m_EventData->GetFirstStage_YE_StripNbr(i), Energy);
+ }
+ }
+ }
+ // First Stage Time
+ unsigned int mysize = m_EventData->GetFirstStageMultXTime();
+ for (UShort_t i = 0; i < mysize; ++i) {
+ Double_t Time= Cal->ApplyCalibration("PISTA/TIME"+NPL::itoa(m_EventData->GetFirstStage_XT_DetectorNbr(i)),m_EventData->GetFirstStage_XT_Time(i));
+ m_PreTreatedData->SetFirstStageXT(m_EventData->GetFirstStage_XT_DetectorNbr(i), m_EventData->GetFirstStage_XT_StripNbr(i), Time);
+ }
+
+ //////
+ // Second Stage Energy
+ sizeFront = m_EventData->GetSecondStageMultXEnergy();
+ for (UShort_t i = 0; i < sizeFront ; ++i) {
+ if (m_EventData->GetSecondStage_XE_Energy(i) > m_E_RAW_Threshold) {
+ Double_t Energy = m_EventData->GetSecondStage_XE_Energy(i);
+ //Double_t Energy = Cal->ApplyCalibration("PISTA/ENERGY"+NPL::itoa(m_EventData->GetSecondStage_XE_DetectorNbr(i)),m_EventData->GetSecondStage_XE_Energy(i));
+ if (Energy > m_E_Threshold) {
+ m_PreTreatedData->SetSecondStageXE(m_EventData->GetSecondStage_XE_DetectorNbr(i), m_EventData->GetSecondStage_XE_StripNbr(i), Energy);
+ }
+ }
+ }
+ sizeBack = m_EventData->GetSecondStageMultXEnergy();
+ for (UShort_t i = 0; i < sizeBack ; ++i) {
+ if (m_EventData->GetSecondStage_YE_Energy(i) > m_E_RAW_Threshold) {
+ Double_t Energy = m_EventData->GetSecondStage_YE_Energy(i);
+ //Double_t Energy = Cal->ApplyCalibration("PISTA/ENERGY"+NPL::itoa(m_EventData->GetSecondStage_YE_DetectorNbr(i)),m_EventData->GetSecondStage_YE_Energy(i));
+ if (Energy > m_E_Threshold) {
+ m_PreTreatedData->SetSecondStageYE(m_EventData->GetSecondStage_YE_DetectorNbr(i), m_EventData->GetSecondStage_YE_StripNbr(i), Energy);
+ }
+ }
+ }
+
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::ReadAnalysisConfig() {
+ bool ReadingStatus = false;
+
+ // path to file
+ string FileName = "./configs/ConfigPISTA.dat";
+
+ // open analysis config file
+ ifstream AnalysisConfigFile;
+ AnalysisConfigFile.open(FileName.c_str());
+
+ if (!AnalysisConfigFile.is_open()) {
+ cout << " No ConfigPISTA.dat found: Default parameter loaded for Analayis " << FileName << endl;
+ return;
+ }
+ cout << " Loading user parameter for Analysis from ConfigPISTA.dat " << endl;
+
+ // Save it in a TAsciiFile
+ TAsciiFile* asciiConfig = RootOutput::getInstance()->GetAsciiFileAnalysisConfig();
+ asciiConfig->AppendLine("%%% ConfigPISTA.dat %%%");
+ asciiConfig->Append(FileName.c_str());
+ asciiConfig->AppendLine("");
+ // read analysis config file
+ string LineBuffer,DataBuffer,whatToDo;
+ while (!AnalysisConfigFile.eof()) {
+ // Pick-up next line
+ getline(AnalysisConfigFile, LineBuffer);
+
+ // search for "header"
+ string name = "ConfigPISTA";
+ if (LineBuffer.compare(0, name.length(), name) == 0)
+ ReadingStatus = true;
+
+ // loop on tokens and data
+ while (ReadingStatus ) {
+ whatToDo="";
+ AnalysisConfigFile >> whatToDo;
+
+ // Search for comment symbol (%)
+ if (whatToDo.compare(0, 1, "%") == 0) {
+ AnalysisConfigFile.ignore(numeric_limits<streamsize>::max(), '\n' );
+ }
+
+ else if (whatToDo=="E_RAW_THRESHOLD") {
+ AnalysisConfigFile >> DataBuffer;
+ m_E_RAW_Threshold = atof(DataBuffer.c_str());
+ cout << whatToDo << " " << m_E_RAW_Threshold << endl;
+ }
+
+ else if (whatToDo=="E_THRESHOLD") {
+ AnalysisConfigFile >> DataBuffer;
+ m_E_Threshold = atof(DataBuffer.c_str());
+ cout << whatToDo << " " << m_E_Threshold << endl;
+ }
+
+ else {
+ ReadingStatus = false;
+ }
+ }
+ }
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::Clear() {
+ EventMultiplicity = 0;
+
+ // Position Information
+ PosX.clear();
+ PosY.clear();
+ PosZ.clear();
+
+ // DSSD
+ DetectorNumber.clear();
+ E.clear();
+ StripX.clear();
+ StripY.clear();
+ Time.clear();
+ DE.clear();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::ReadConfiguration(NPL::InputParser parser) {
+ vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("PISTA");
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << "//// " << blocks.size() << " detectors found " << endl;
+
+ vector<string> cart = {"POS"};
+ vector<string> sphe = {"R","Theta","Phi"};
+
+ for(unsigned int i = 0 ; i < blocks.size() ; i++){
+ if(blocks[i]->HasTokenList(cart)){
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// PISTA " << i+1 << endl;
+
+ TVector3 Pos = blocks[i]->GetTVector3("POS","mm");
+
+ AddDetector(Pos);
+ }
+ else if(blocks[i]->HasTokenList(sphe)){
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// PISTA " << i+1 << endl;
+ double R = blocks[i]->GetDouble("R","mm");
+ double Theta = blocks[i]->GetDouble("Theta","deg");
+ double Phi = blocks[i]->GetDouble("Phi","deg");
+
+ AddDetector(R, Theta, Phi);
+ }
+ else{
+ cout << "ERROR: check your input file formatting " << endl;
+ exit(1);
+ }
+ }
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::InitSpectra() {
+ m_Spectra = new TPISTASpectra(m_NumberOfDetectors);
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::FillSpectra() {
+ m_Spectra -> FillRawSpectra(m_EventData);
+ m_Spectra -> FillPreTreatedSpectra(m_PreTreatedData);
+ m_Spectra -> FillPhysicsSpectra(m_EventPhysics);
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::CheckSpectra() {
+ m_Spectra->CheckSpectra();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::ClearSpectra() {
+ // To be done
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+map< string , TH1*> TPISTAPhysics::GetSpectra() {
+ if(m_Spectra)
+ return m_Spectra->GetMapHisto();
+ else{
+ map< string , TH1*> empty;
+ return empty;
+ }
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::WriteSpectra() {
+ m_Spectra->WriteSpectra();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::AddParameterToCalibrationManager() {
+ CalibrationManager* Cal = CalibrationManager::getInstance();
+ for (int i = 0; i < m_NumberOfDetectors; ++i) {
+ Cal->AddParameter("PISTA", "D"+ NPL::itoa(i+1)+"_ENERGY","PISTA_D"+ NPL::itoa(i+1)+"_ENERGY");
+ Cal->AddParameter("PISTA", "D"+ NPL::itoa(i+1)+"_TIME","PISTA_D"+ NPL::itoa(i+1)+"_TIME");
+ }
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::InitializeRootInputRaw() {
+ TChain* inputChain = RootInput::getInstance()->GetChain();
+ inputChain->SetBranchStatus("PISTA", true );
+ inputChain->SetBranchAddress("PISTA", &m_EventData );
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::InitializeRootInputPhysics() {
+ TChain* inputChain = RootInput::getInstance()->GetChain();
+ inputChain->SetBranchAddress("PISTA", &m_EventPhysics);
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TPISTAPhysics::InitializeRootOutput() {
+ TTree* outputTree = RootOutput::getInstance()->GetTree();
+ outputTree->Branch("PISTA", "TPISTAPhysics", &m_EventPhysics);
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+// Construct Method to be pass to the DetectorFactory //
+////////////////////////////////////////////////////////////////////////////////
+NPL::VDetector* TPISTAPhysics::Construct() {
+ return (NPL::VDetector*) new TPISTAPhysics();
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+// Registering the construct method to the factory //
+////////////////////////////////////////////////////////////////////////////////
+extern "C"{
+ class proxy_PISTA{
+ public:
+ proxy_PISTA(){
+ NPL::DetectorFactory::getInstance()->AddToken("PISTA","PISTA");
+ NPL::DetectorFactory::getInstance()->AddDetector("PISTA",TPISTAPhysics::Construct);
+ }
+ };
+
+ proxy_PISTA p_PISTA;
+}
+
diff --git a/NPLib/Detectors/PISTA/TPISTAPhysics.h b/NPLib/Detectors/PISTA/TPISTAPhysics.h
new file mode 100644
index 0000000000000000000000000000000000000000..e2c2071a49957200e2c02753f067c4b08ebbbe11
--- /dev/null
+++ b/NPLib/Detectors/PISTA/TPISTAPhysics.h
@@ -0,0 +1,216 @@
+#ifndef TPISTAPHYSICS_H
+#define TPISTAPHYSICS_H
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold PISTA Treated data *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+// C++ headers
+#include <vector>
+#include <map>
+#include <string>
+using namespace std;
+
+// ROOT headers
+#include "TObject.h"
+#include "TH1.h"
+#include "TVector3.h"
+// NPTool headers
+#include "TPISTAData.h"
+#include "TPISTASpectra.h"
+#include "NPCalibrationManager.h"
+#include "NPVDetector.h"
+#include "NPInputParser.h"
+// forward declaration
+class TPISTASpectra;
+
+
+
+class TPISTAPhysics : public TObject, public NPL::VDetector {
+ //////////////////////////////////////////////////////////////
+ // constructor and destructor
+ public:
+ TPISTAPhysics();
+ ~TPISTAPhysics() {};
+
+
+ //////////////////////////////////////////////////////////////
+ // Inherited from TObject and overriden to avoid warnings
+ public:
+ void Clear();
+ void Clear(const Option_t*) {};
+
+ public:
+ vector<TVector2> Match_X_Y();
+ int CheckEvent();
+
+ //////////////////////////////////////////////////////////////
+ // data obtained after BuildPhysicalEvent() and stored in
+ // output ROOT file
+ public:
+ Int_t EventMultiplicity;
+ vector<int> DetectorNumber;
+ vector<double> E;
+ vector<double> DE;
+ vector<int> StripX;
+ vector<int> StripY;
+ vector<double> Time;
+
+ vector<double> PosX;
+ vector<double> PosY;
+ vector<double> PosZ;
+
+
+ //////////////////////////////////////////////////////////////
+ // methods inherited from the VDetector ABC class
+ public:
+ // read stream from ConfigFile to pick-up detector parameters
+ void ReadConfiguration(NPL::InputParser);
+
+ /// A usefull method to bundle all operation to add a detector
+ void AddDetector(TVector3 POS);
+ void AddDetector(double R, double Theta, double Phi);
+
+ // add parameters to the CalibrationManger
+ void AddParameterToCalibrationManager();
+
+ // method called event by event, aiming at extracting the
+ // physical information from detector
+ void BuildPhysicalEvent();
+
+ // same as BuildPhysicalEvent() method but with a simpler
+ // treatment
+ void BuildSimplePhysicalEvent();
+
+ // same as above but for online analysis
+ void BuildOnlinePhysicalEvent() {BuildPhysicalEvent();};
+
+ // activate raw data object and branches from input TChain
+ // in this method mother branches (Detector) AND daughter leaves
+ // (fDetector_parameter) have to be activated
+ void InitializeRootInputRaw();
+
+ // activate physics data object and branches from input TChain
+ // in this method mother branches (Detector) AND daughter leaves
+ // (fDetector_parameter) have to be activated
+ void InitializeRootInputPhysics();
+
+ // create branches of output ROOT file
+ void InitializeRootOutput();
+
+ // clear the raw and physical data objects event by event
+ void ClearEventPhysics() {Clear();}
+ void ClearEventData() {m_EventData->Clear();}
+
+ // methods related to the TPISTASpectra class
+ // instantiate the TPISTASpectra class and
+ // declare list of histograms
+ void InitSpectra();
+
+ // fill the spectra
+ void FillSpectra();
+
+ // used for Online mainly, sanity check for histograms and
+ // change their color if issues are found, for example
+ void CheckSpectra();
+
+ // used for Online only, clear all the spectra
+ void ClearSpectra();
+
+ // write spectra to ROOT output file
+ void WriteSpectra();
+
+
+ //////////////////////////////////////////////////////////////
+ // specific methods to PISTA array
+ public:
+ // remove bad channels, calibrate the data and apply thresholds
+ void PreTreat();
+
+ // clear the pre-treated object
+ void ClearPreTreatedData() {m_PreTreatedData->Clear();}
+
+ // read the user configuration file. If no file is found, load standard one
+ void ReadAnalysisConfig();
+
+ // give and external TPISTAData object to TPISTAPhysics.
+ // needed for online analysis for example
+ void SetRawDataPointer(TPISTAData* rawDataPointer) {m_EventData = rawDataPointer;}
+
+ double GetNumberOfTelescope() const {return m_NumberOfDetectors;}
+ int GetEventMultiplicity() const {return EventMultiplicity;}
+
+ double GetStripPositionX(const int N, const int X, const int Y){
+ return m_StripPositionX[N-1][X-1][Y-1];
+ };
+ double GetStripPositionY(const int N, const int X, const int Y){
+ return m_StripPositionY[N-1][X-1][Y-1];
+ };
+ double GetStripPositionZ(const int N, const int X, const int Y){
+ return m_StripPositionZ[N-1][X-1][Y-1];
+ };
+
+
+ TVector3 GetPositionOfInteraction(const int i);
+ TVector3 GetDetectorNormal(const int i);
+
+ // objects are not written in the TTree
+ private:
+ TPISTAData* m_EventData; //!
+ TPISTAData* m_PreTreatedData; //!
+ TPISTAPhysics* m_EventPhysics; //!
+
+ // getters for raw and pre-treated data object
+ public:
+ TPISTAData* GetRawData() const {return m_EventData;}
+ TPISTAData* GetPreTreatedData() const {return m_PreTreatedData;}
+
+ // parameters used in the analysis
+ private:
+ int m_NumberOfDetectors; //!
+
+ vector<vector<vector<double>>> m_StripPositionX; //!
+ vector<vector<vector<double>>> m_StripPositionY; //!
+ vector<vector<vector<double>>> m_StripPositionZ; //!
+
+ // thresholds
+ double m_E_RAW_Threshold; //!
+ double m_E_Threshold; //!
+
+ private:
+ unsigned int m_MaximumStripMultiplicityAllowed;//
+ double m_StripEnergyMatching;//
+
+
+ // spectra class
+ private:
+ TPISTASpectra* m_Spectra; // !
+
+ // spectra getter
+ public:
+ map<string, TH1*> GetSpectra();
+
+ // Static constructor to be passed to the Detector Factory
+ public:
+ static NPL::VDetector* Construct();
+
+ ClassDef(TPISTAPhysics,1) // PISTAPhysics structure
+};
+#endif
diff --git a/NPLib/Detectors/PISTA/TPISTASpectra.cxx b/NPLib/Detectors/PISTA/TPISTASpectra.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..c40fe4d4fe86b5f494c4255f199987b878c8e0e1
--- /dev/null
+++ b/NPLib/Detectors/PISTA/TPISTASpectra.cxx
@@ -0,0 +1,174 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold PISTA Spectra *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+// class header
+#include "TPISTASpectra.h"
+
+// STL
+#include <iostream>
+#include <string>
+using namespace std;
+
+// NPTool header
+#include "NPOptionManager.h"
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+TPISTASpectra::TPISTASpectra()
+ : fNumberOfDetectors(0) {
+ SetName("PISTA");
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+TPISTASpectra::TPISTASpectra(unsigned int NumberOfDetectors) {
+ if(NPOptionManager::getInstance()->GetVerboseLevel()>0)
+ cout << "************************************************" << endl
+ << "TPISTASpectra : Initalizing control spectra for "
+ << NumberOfDetectors << " Detectors" << endl
+ << "************************************************" << endl ;
+ SetName("PISTA");
+ fNumberOfDetectors = NumberOfDetectors;
+
+ InitRawSpectra();
+ InitPreTreatedSpectra();
+ InitPhysicsSpectra();
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+TPISTASpectra::~TPISTASpectra() {
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TPISTASpectra::InitRawSpectra() {
+ static string name;
+ for (unsigned int i = 0; i < fNumberOfDetectors; i++) { // loop on number of detectors
+ // Energy
+ name = "PISTA"+NPL::itoa(i+1)+"_ENERGY_RAW";
+ AddHisto1D(name, name, 4096, 0, 16384, "PISTA/RAW");
+ // Time
+ name = "PISTA"+NPL::itoa(i+1)+"_TIME_RAW";
+ AddHisto1D(name, name, 4096, 0, 16384, "PISTA/RAW");
+ } // end loop on number of detectors
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TPISTASpectra::InitPreTreatedSpectra() {
+ static string name;
+ for (unsigned int i = 0; i < fNumberOfDetectors; i++) { // loop on number of detectors
+ // Energy
+ name = "PISTA"+NPL::itoa(i+1)+"_ENERGY_CAL";
+ AddHisto1D(name, name, 500, 0, 25, "PISTA/CAL");
+ // Time
+ name = "PISTA"+NPL::itoa(i+1)+"_TIME_CAL";
+ AddHisto1D(name, name, 500, 0, 25, "PISTA/CAL");
+
+
+ } // end loop on number of detectors
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TPISTASpectra::InitPhysicsSpectra() {
+ static string name;
+ // Kinematic Plot
+ name = "PISTA_ENERGY_TIME";
+ AddHisto2D(name, name, 500, 0, 500, 500, 0, 50, "PISTA/PHY");
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TPISTASpectra::FillRawSpectra(TPISTAData* RawData) {
+ static string name;
+ static string family;
+
+ // Energy
+ unsigned int sizeE = RawData->GetFirstStageMultXEnergy();
+ for (unsigned int i = 0; i < sizeE; i++) {
+ name = "PISTA"+NPL::itoa(RawData->GetFirstStage_XE_DetectorNbr(i))+"_ENERGY_RAW";
+ family = "PISTA/RAW";
+
+ FillSpectra(family,name,RawData->GetFirstStage_XE_Energy(i));
+ }
+
+ // Time
+ unsigned int sizeT = RawData->GetFirstStageMultXTime();
+ for (unsigned int i = 0; i < sizeT; i++) {
+ name = "PISTA"+NPL::itoa(RawData->GetFirstStage_XT_DetectorNbr(i))+"_TIME_RAW";
+ family = "PISTA/RAW";
+
+ FillSpectra(family,name,RawData->GetFirstStage_XT_Time(i));
+ }
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TPISTASpectra::FillPreTreatedSpectra(TPISTAData* PreTreatedData) {
+ static string name;
+ static string family;
+
+ // Energy
+ unsigned int sizeE = PreTreatedData->GetFirstStageMultXEnergy();
+ for (unsigned int i = 0; i < sizeE; i++) {
+ name = "PISTA"+NPL::itoa(PreTreatedData->GetFirstStage_XE_DetectorNbr(i))+"_ENERGY_CAL";
+ family = "PISTA/CAL";
+
+ FillSpectra(family,name,PreTreatedData->GetFirstStage_XE_Energy(i));
+ }
+
+ // Time
+ unsigned int sizeT = PreTreatedData->GetFirstStageMultXTime();
+ for (unsigned int i = 0; i < sizeT; i++) {
+ name = "PISTA"+NPL::itoa(PreTreatedData->GetFirstStage_XT_DetectorNbr(i))+"_TIME_CAL";
+ family = "PISTA/CAL";
+
+ FillSpectra(family,name,PreTreatedData->GetFirstStage_XT_Time(i));
+ }
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TPISTASpectra::FillPhysicsSpectra(TPISTAPhysics* Physics) {
+ static string name;
+ static string family;
+ family= "PISTA/PHY";
+
+ // Energy vs time
+ unsigned int sizeE = Physics->E.size();
+ for(unsigned int i = 0 ; i < sizeE ; i++){
+ name = "PISTA_ENERGY_TIME";
+ FillSpectra(family,name,Physics->E[i],Physics->Time[i]);
+ }
+}
+
diff --git a/NPLib/Detectors/PISTA/TPISTASpectra.h b/NPLib/Detectors/PISTA/TPISTASpectra.h
new file mode 100644
index 0000000000000000000000000000000000000000..97239755cd2efb3e9216c0caab255aed89d004ad
--- /dev/null
+++ b/NPLib/Detectors/PISTA/TPISTASpectra.h
@@ -0,0 +1,62 @@
+#ifndef TPISTASPECTRA_H
+#define TPISTASPECTRA_H
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold PISTA Spectra *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+// NPLib headers
+#include "NPVSpectra.h"
+#include "TPISTAData.h"
+#include "TPISTAPhysics.h"
+
+// Forward Declaration
+class TPISTAPhysics;
+
+
+class TPISTASpectra : public VSpectra {
+ //////////////////////////////////////////////////////////////
+ // constructor and destructor
+ public:
+ TPISTASpectra();
+ TPISTASpectra(unsigned int NumberOfDetectors);
+ ~TPISTASpectra();
+
+ //////////////////////////////////////////////////////////////
+ // Initialization methods
+ private:
+ void InitRawSpectra();
+ void InitPreTreatedSpectra();
+ void InitPhysicsSpectra();
+
+ //////////////////////////////////////////////////////////////
+ // Filling methods
+ public:
+ void FillRawSpectra(TPISTAData*);
+ void FillPreTreatedSpectra(TPISTAData*);
+ void FillPhysicsSpectra(TPISTAPhysics*);
+
+ //////////////////////////////////////////////////////////////
+ // Detector parameters
+ private:
+ unsigned int fNumberOfDetectors;
+};
+
+#endif
diff --git a/NPLib/Detectors/Scone/TSconePhysics.cxx b/NPLib/Detectors/Scone/TSconePhysics.cxx
index b30e95cf0bf8011fc6a016477a99aac76dcbce1c..54031fe4c5e4cb5c4729a0b0e274b8de92f9d6dd 100644
--- a/NPLib/Detectors/Scone/TSconePhysics.cxx
+++ b/NPLib/Detectors/Scone/TSconePhysics.cxx
@@ -87,11 +87,13 @@ void TSconePhysics::BuildPhysicalEvent() {
// match energy and time together
unsigned int mysizeE = m_PreTreatedData->GetMultEnergy();
unsigned int mysizeT = m_PreTreatedData->GetMultTime();
- for (UShort_t e = 0; e < mysizeE ; e++) {
+ //for (UShort_t e = 0; e < mysizeE ; e++) {
+ if (mysizeE == mysizeT) {
for (UShort_t t = 0; t < mysizeT ; t++) {
- if (m_PreTreatedData->GetE_DetectorNbr(e) == m_PreTreatedData->GetT_DetectorNbr(t)) {
- int det = m_PreTreatedData->GetE_DetectorNbr(e);
- Energy_det[det-1] += m_PreTreatedData->Get_Energy(e);
+ if (m_PreTreatedData->GetE_DetectorNbr(t) == m_PreTreatedData->GetT_DetectorNbr(t)) {
+ int det = m_PreTreatedData->GetE_DetectorNbr(t);
+ Energy_det[det-1] += m_PreTreatedData->Get_Energy(t);
+ //cout << det << " " << Energy_det[det-1] << " " << m_PreTreatedData->GetE_PlasticNbr(t) << endl;
if(m_PreTreatedData->Get_Time(t)>Time_det[det-1]){
Time_det[det-1] = m_PreTreatedData->Get_Time(t);
}
diff --git a/NPLib/Physics/NPReaction.cxx b/NPLib/Physics/NPReaction.cxx
index ba6e8803ba0f5c3638a47d78c111eaf841e4d9c2..566194f56458581b04e8dd22f16548e31bec6c4c 100644
--- a/NPLib/Physics/NPReaction.cxx
+++ b/NPLib/Physics/NPReaction.cxx
@@ -87,7 +87,6 @@ ClassImp(Reaction)
fshoot3=true;
fshoot4=true;
fUseExInGeant4=true;
- RandGen=gRandom;
fLabCrossSection=false; // flag if the provided cross-section is in the lab or not
@@ -156,7 +155,6 @@ Reaction::Reaction(string reaction){
fshoot3=true;
fshoot4=true;
- RandGen=gRandom;
fLabCrossSection=false;
@@ -357,7 +355,7 @@ double Reaction::EnergyLabFromThetaLab(double ThetaLab){
if(B>D) {
ThetaLabMax = asin(sqrt(D/B));
- if(RandGen->Rndm()<0.5) sign=-1;
+ if(gRandom->Rndm()<0.5) sign=-1;
}
if(ThetaLab>ThetaLabMax) return -1;
diff --git a/NPSimulation/Core/MaterialManager.cc b/NPSimulation/Core/MaterialManager.cc
index 5ad7f38596a949a074331483574803cc30e3060e..e66bc6168a3e9aa0269523fcf3bf3732af3b4bd2 100644
--- a/NPSimulation/Core/MaterialManager.cc
+++ b/NPSimulation/Core/MaterialManager.cc
@@ -135,6 +135,19 @@ G4Material* MaterialManager::GetMaterialFromLibrary(string Name,
return material;
}
+ else if (Name == "Rogers4003C") {
+ if (!density)
+ density = 1.79 * g / cm3;
+ G4Material* material = new G4Material("NPS_" + Name, density, 4);
+ material->AddElement(GetElementFromLibrary("H"), 2);
+ material->AddElement(GetElementFromLibrary("C"), 50);
+ material->AddElement(GetElementFromLibrary("O"), 38);
+ material->AddElement(GetElementFromLibrary("Si"), 10);
+ m_Material[Name] = material;
+ return material;
+ }
+
+
else if (Name == "Mylar") {
if (!density)
density = 1.397 * g / cm3;
@@ -1192,6 +1205,7 @@ void MaterialManager::WriteCrossSectionTable(G4ParticleDefinition* Particle, G4d
G4String path1;
G4String path2;
G4String path3;
+ G4String path4;
G4String ParticleName = Particle->GetParticleName();
G4String MaterialName = it->second->GetName();
G4String ElementName = it->second->GetElement(i)->GetName();
@@ -1199,13 +1213,16 @@ void MaterialManager::WriteCrossSectionTable(G4ParticleDefinition* Particle, G4d
path1 = GlobalPath + "/Inputs/CrossSection/" + "G4XS_elastic_" + ParticleName + "_" + ElementName + ".dat";
path2 = GlobalPath + "/Inputs/CrossSection/" + "G4XS_inelastic_" + ParticleName + "_" + ElementName + ".dat";
path3 = GlobalPath + "/Inputs/CrossSection/" + "G4XS_capture_" + ParticleName + "_" + ElementName + ".dat";
+ path4 = GlobalPath + "/Inputs/CrossSection/" + "G4XS_fission_" + ParticleName + "_" + ElementName + ".dat";
ofstream ofile_elastic;
ofstream ofile_inelastic;
ofstream ofile_capture;
+ ofstream ofile_fission;
ofile_elastic.open(path1);
ofile_inelastic.open(path2);
ofile_capture.open(path3);
+ ofile_fission.open(path4);
//std::cout << path << std::endl;
double xs;
double step_keV = 1*keV;
@@ -1230,12 +1247,16 @@ void MaterialManager::WriteCrossSectionTable(G4ParticleDefinition* Particle, G4d
// Capture Cross Section
xs = store->GetCaptureCrossSectionPerAtom(Particle, E, it->second->GetElement(i), it->second);
ofile_capture << E/MeV << " " << xs/barn << G4endl;
-
+
+ // Fission Cross Section
+ xs = store->GetFissionCrossSectionPerAtom(Particle, E, it->second->GetElement(i), it->second);
+ ofile_fission << E/MeV << " " << xs/barn << G4endl;
}
ofile_elastic.close();
ofile_inelastic.close();
ofile_capture.close();
+ ofile_fission.close();
}
}
}
diff --git a/NPSimulation/Detectors/ChiNu/ChiNu.cc b/NPSimulation/Detectors/ChiNu/ChiNu.cc
index 7d747ef2a82a010a5e40e0a58c9323db06d22523..1f948124ecc42431b681979acc051ec1cba2b5fe 100644
--- a/NPSimulation/Detectors/ChiNu/ChiNu.cc
+++ b/NPSimulation/Detectors/ChiNu/ChiNu.cc
@@ -444,10 +444,10 @@ void ChiNu::ConstructDetector(G4LogicalVolume* world){
BuildDetector()->MakeImprint(world,Det_pos,Rot,i);
}
- /*G4RotationMatrix* Rot_FC = new G4RotationMatrix(0,0,0);
+ G4RotationMatrix* Rot_FC = new G4RotationMatrix(0,0,0);
G4ThreeVector Pos_FC = G4ThreeVector(0,0,0) ;
BuildFissionChamber()->MakeImprint(world,Pos_FC,Rot_FC,0);
-*/
+
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
diff --git a/NPSimulation/Detectors/PISTA/CMakeLists.txt b/NPSimulation/Detectors/PISTA/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..475afe7666c80dcdfd88a6916cc4cc8e12a4dfaf
--- /dev/null
+++ b/NPSimulation/Detectors/PISTA/CMakeLists.txt
@@ -0,0 +1,2 @@
+add_library(NPSPISTA SHARED PISTA.cc)
+target_link_libraries(NPSPISTA NPSCore ${ROOT_LIBRARIES} ${Geant4_LIBRARIES} ${NPLib_LIBRARIES} -lNPPISTA)
diff --git a/NPSimulation/Detectors/PISTA/PISTA.cc b/NPSimulation/Detectors/PISTA/PISTA.cc
new file mode 100644
index 0000000000000000000000000000000000000000..0a02fb280cd0201f8d7fa0f46c2ee74e747111bf
--- /dev/null
+++ b/NPSimulation/Detectors/PISTA/PISTA.cc
@@ -0,0 +1,387 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class describe PISTA simulation *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ *****************************************************************************/
+
+// C++ headers
+#include <sstream>
+#include <cmath>
+#include <limits>
+//G4 Geometry object
+#include "G4Tubs.hh"
+#include "G4Trap.hh"
+#include "G4Box.hh"
+
+//G4 sensitive
+#include "G4SDManager.hh"
+#include "G4MultiFunctionalDetector.hh"
+
+//G4 various object
+#include "G4Material.hh"
+#include "G4Transform3D.hh"
+#include "G4PVPlacement.hh"
+#include "G4VisAttributes.hh"
+#include "G4Colour.hh"
+
+// NPTool header
+#include "PISTA.hh"
+#include "DSSDScorers.hh"
+#include "InteractionScorers.hh"
+#include "RootOutput.h"
+#include "MaterialManager.hh"
+#include "NPSDetectorFactory.hh"
+#include "NPOptionManager.h"
+#include "NPSHitsMap.hh"
+// CLHEP header
+#include "CLHEP/Random/RandGauss.h"
+
+using namespace std;
+using namespace CLHEP;
+
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+namespace PISTA_NS{
+ // Energy and time Resolution
+ const double EnergyThreshold = 0.1*MeV;
+ const double ResoTime = 0.2*ns ;
+ const double ResoEnergy = 0.015*MeV ;
+
+ // Trapezoid dimension
+ const double TrapezoidBaseLarge = 95*mm;
+ const double TrapezoidBaseSmall = 45*mm;
+ const double TrapezoidHeight = 118*mm;
+ const double TrapezoidLength = 1*cm;
+ const double FirstStageThickness = 100*um;
+ const double SecondStageThickness = 1*mm;
+ const double DistanceBetweenSi = 7*mm;
+ //const double FirstStageNbrOfStrips = 128;
+ //const double SecondStageNbrOfStrips = 16;
+}
+using namespace PISTA_NS;
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+// PISTA Specific Method
+PISTA::PISTA(){
+ m_Event = new TPISTAData() ;
+ m_FirstStageScorer = 0;
+ m_SecondStageScorer = 0;
+ m_TrapezoidDetector = 0;
+}
+
+PISTA::~PISTA(){
+}
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+void PISTA::AddDetector(G4ThreeVector POS){
+ // Convert the POS value to R theta Phi as Spherical coordinate is easier in G4
+ m_R.push_back(POS.mag());
+ m_Theta.push_back(POS.theta());
+ m_Phi.push_back(POS.phi());
+}
+
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+void PISTA::AddDetector(double R, double Theta, double Phi){
+ m_R.push_back(R);
+ m_Theta.push_back(Theta);
+ m_Phi.push_back(Phi);
+}
+
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+G4LogicalVolume* PISTA::BuildTrapezoidDetector(){
+ if(!m_TrapezoidDetector){
+ // Definittion of the volume containing the sensitive detectors
+ G4Trap* solidTrapezoid = new G4Trap("PISTA",
+ TrapezoidLength*0.5, 0*deg, 0*deg,
+ TrapezoidHeight*0.5, TrapezoidBaseLarge*0.5, TrapezoidBaseSmall*0.5, 0*deg,
+ TrapezoidHeight*0.5, TrapezoidBaseLarge*0.5, TrapezoidBaseSmall*0.5, 0*deg);
+ G4LogicalVolume* logicTrapezoid = new G4LogicalVolume(solidTrapezoid,
+ MaterialManager::getInstance()->GetMaterialFromLibrary("Vacuum"),
+ "PISTA",
+ 0,0,0);
+
+ G4VisAttributes* TrapezoidVisAtt = new G4VisAttributes(G4Colour(0.90, 0.90, 0.90));
+ TrapezoidVisAtt->SetForceWireframe(true);
+ logicTrapezoid->SetVisAttributes(TrapezoidVisAtt);
+
+ // First stage silicon detector
+ G4ThreeVector positionFirstStage = G4ThreeVector(0,0,-4*mm);
+
+ G4Trap* solidFirstStage = new G4Trap("solidFirstSatge",
+ FirstStageThickness*0.5, 0*deg, 0*deg,
+ TrapezoidHeight*0.5, TrapezoidBaseLarge*0.5, TrapezoidBaseSmall*0.5, 0*deg,
+ TrapezoidHeight*0.5, TrapezoidBaseLarge*0.5, TrapezoidBaseSmall*0.5, 0*deg);
+ G4LogicalVolume* logicFirstStage = new G4LogicalVolume(solidFirstStage,
+ MaterialManager::getInstance()->GetMaterialFromLibrary("Si"),
+ "logicFirstStage",
+ 0,0,0);
+ new G4PVPlacement(0,
+ positionFirstStage,
+ logicFirstStage,
+ "PISTA_FirstStage",
+ logicTrapezoid,
+ false,
+ 0);
+ // Set First Stage sensitive
+ logicFirstStage->SetSensitiveDetector(m_FirstStageScorer);
+
+ // Visualisation of First Stage strips
+ G4VisAttributes* FirstStageVisAtt = new G4VisAttributes(G4Colour(0.3,0.3,0.3));
+ logicFirstStage->SetVisAttributes(FirstStageVisAtt);
+
+ //////
+ // Second stage silicon detector
+ G4ThreeVector positionSecondStage = G4ThreeVector(0,0,-0.5*TrapezoidLength+DistanceBetweenSi);
+
+ G4Trap* solidSecondStage = new G4Trap("solidSecondSatge",
+ SecondStageThickness*0.5, 0*deg, 0*deg,
+ TrapezoidHeight*0.5, TrapezoidBaseLarge*0.5, TrapezoidBaseSmall*0.5, 0*deg,
+ TrapezoidHeight*0.5, TrapezoidBaseLarge*0.5, TrapezoidBaseSmall*0.5, 0*deg);
+ G4LogicalVolume* logicSecondStage = new G4LogicalVolume(solidSecondStage,
+ MaterialManager::getInstance()->GetMaterialFromLibrary("Si"),
+ "logicSecondStage",
+ 0,0,0);
+ new G4PVPlacement(0,
+ positionSecondStage,
+ logicSecondStage,
+ "PISTA_SecondStage",
+ logicTrapezoid,
+ false,
+ 0);
+ // Set Second Stage sensitive
+ logicSecondStage->SetSensitiveDetector(m_SecondStageScorer);
+
+ // Visualisation of Second Stage strips
+ G4VisAttributes* SecondStageVisAtt = new G4VisAttributes(G4Colour(0.4,0.5,0.5));
+ logicSecondStage->SetVisAttributes(SecondStageVisAtt);
+
+
+ m_TrapezoidDetector = logicTrapezoid;
+ }
+ return m_TrapezoidDetector;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+// Virtual Method of NPS::VDetector class
+
+// Read stream at Configfile to pick-up parameters of detector (Position,...)
+// Called in DetecorConstruction::ReadDetextorConfiguration Method
+void PISTA::ReadConfiguration(NPL::InputParser parser){
+ vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("PISTA");
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << "//// " << blocks.size() << " detectors found " << endl;
+
+ vector<string> cart = {"POS"};
+ vector<string> sphe = {"R","Theta","Phi"};
+
+ for(unsigned int i = 0 ; i < blocks.size() ; i++){
+ if(blocks[i]->HasTokenList(cart)){
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// PISTA " << i+1 << endl;
+
+ G4ThreeVector Pos = NPS::ConvertVector(blocks[i]->GetTVector3("POS","mm"));
+ AddDetector(Pos);
+ }
+ else if(blocks[i]->HasTokenList(sphe)){
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// PISTA " << i+1 << endl;
+ double R = blocks[i]->GetDouble("R","mm");
+ double Theta = blocks[i]->GetDouble("Theta","deg");
+ double Phi = blocks[i]->GetDouble("Phi","deg");
+ AddDetector(R,Theta,Phi);
+ }
+ else{
+ cout << "ERROR: check your input file formatting " << endl;
+ exit(1);
+ }
+ }
+}
+
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
+// Construct detector and inialise sensitive part.
+// Called After DetecorConstruction::AddDetector Method
+void PISTA::ConstructDetector(G4LogicalVolume* world){
+ for (unsigned short i = 0 ; i < m_R.size() ; i++) {
+
+ G4double wX = m_R[i] * sin(m_Theta[i] ) * cos(m_Phi[i] ) ;
+ G4double wY = m_R[i] * sin(m_Theta[i] ) * sin(m_Phi[i] ) ;
+ G4double wZ = m_R[i] * cos(m_Theta[i] ) ;
+ G4ThreeVector Det_pos = G4ThreeVector(wX, wY, wZ) ;
+ // So the face of the detector is at R instead of the middle
+ Det_pos+=Det_pos.unit()*PISTA_NS::TrapezoidLength*0.5;
+ // Building Detector reference frame
+ G4double ii = cos(m_Theta[i]) * cos(m_Phi[i]);
+ G4double jj = cos(m_Theta[i]) * sin(m_Phi[i]);
+ G4double kk = -sin(m_Theta[i]);
+ G4ThreeVector Y(ii,jj,kk);
+ G4ThreeVector w = Det_pos.unit();
+ G4ThreeVector u = w.cross(Y);
+ G4ThreeVector v = w.cross(u);
+ v = v.unit();
+ u = u.unit();
+
+ G4RotationMatrix* Rot = new G4RotationMatrix(u,v,w);
+
+ new G4PVPlacement(G4Transform3D(*Rot,Det_pos),
+ BuildTrapezoidDetector(),
+ "PISTA",world,false,i+1);
+
+ }
+}
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+// Add Detector branch to the EventTree.
+// Called After DetecorConstruction::AddDetector Method
+void PISTA::InitializeRootOutput(){
+ RootOutput *pAnalysis = RootOutput::getInstance();
+ TTree *pTree = pAnalysis->GetTree();
+ if(!pTree->FindBranch("PISTA")){
+ pTree->Branch("PISTA", "TPISTAData", &m_Event) ;
+ }
+ pTree->SetBranchAddress("PISTA", &m_Event) ;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+// Read sensitive part and fill the Root tree.
+// Called at in the EventAction::EndOfEventAvtion
+void PISTA::ReadSensitive(const G4Event* ){
+ m_Event->Clear();
+
+ ///////////
+ // First Stage scorer
+ DSSDScorers::PS_Rectangle* FirstStageScorer= (DSSDScorers::PS_Rectangle*) m_FirstStageScorer->GetPrimitive(0);
+
+ unsigned int sizeFront = FirstStageScorer->GetLengthMult();
+ for(unsigned int i = 0 ; i < sizeFront ; i++){
+ double Energy = RandGauss::shoot(FirstStageScorer->GetEnergyLength(i), ResoEnergy);
+ if(Energy>EnergyThreshold){
+ double Time = RandGauss::shoot(FirstStageScorer->GetTimeLength(i), ResoTime);
+ int DetNbr = FirstStageScorer->GetDetectorLength(i);
+ int StripFront = FirstStageScorer->GetStripLength(i);
+ m_Event->SetFirstStageXE(DetNbr, StripFront, Energy);
+ m_Event->SetFirstStageXT(DetNbr, StripFront, Time);
+ }
+ }
+ unsigned int sizeBack = FirstStageScorer->GetWidthMult();
+ for(unsigned int i = 0 ; i < sizeBack ; i++){
+ double Energy = RandGauss::shoot(FirstStageScorer->GetEnergyWidth(i), ResoEnergy);
+ if(Energy>EnergyThreshold){
+ double Time = RandGauss::shoot(FirstStageScorer->GetTimeWidth(i), ResoTime);
+ int DetNbr = FirstStageScorer->GetDetectorWidth(i);
+ int StripFront = FirstStageScorer->GetStripWidth(i);
+ m_Event->SetFirstStageYE(DetNbr, StripFront, Energy);
+ m_Event->SetFirstStageYT(DetNbr, StripFront, Time);
+ }
+ }
+ FirstStageScorer->clear();
+
+ ///////////
+ // Second Stage scorer
+ DSSDScorers::PS_Rectangle* SecondStageScorer= (DSSDScorers::PS_Rectangle*) m_SecondStageScorer->GetPrimitive(0);
+
+ unsigned int sizeFrontSecondStage = SecondStageScorer->GetLengthMult();
+ for(unsigned int i = 0 ; i < sizeFrontSecondStage ; i++){
+ double Energy = RandGauss::shoot(SecondStageScorer->GetEnergyLength(i), ResoEnergy);
+ if(Energy>EnergyThreshold){
+ double Time = RandGauss::shoot(SecondStageScorer->GetTimeLength(i), ResoTime);
+ int DetNbr = SecondStageScorer->GetDetectorLength(i);
+ int StripFront = SecondStageScorer->GetStripLength(i);
+ m_Event->SetSecondStageXE(DetNbr, StripFront, Energy);
+ m_Event->SetSecondStageXT(DetNbr, StripFront, Time);
+ }
+ }
+ unsigned int sizeBackSecondStage = SecondStageScorer->GetWidthMult();
+ for(unsigned int i = 0 ; i < sizeBackSecondStage ; i++){
+ double Energy = RandGauss::shoot(SecondStageScorer->GetEnergyWidth(i), ResoEnergy);
+ if(Energy>EnergyThreshold){
+ double Time = RandGauss::shoot(SecondStageScorer->GetTimeWidth(i), ResoTime);
+ int DetNbr = SecondStageScorer->GetDetectorWidth(i);
+ int StripFront = SecondStageScorer->GetStripWidth(i);
+ m_Event->SetSecondStageYE(DetNbr, StripFront, Energy);
+ m_Event->SetSecondStageYT(DetNbr, StripFront, Time);
+ }
+ }
+ SecondStageScorer->clear();
+
+
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+////////////////////////////////////////////////////////////////
+void PISTA::InitializeScorers() {
+ // This check is necessary in case the geometry is reloaded
+ bool already_exist = false;
+ m_FirstStageScorer = CheckScorer("FirstStageScorer",already_exist) ;
+ m_SecondStageScorer = CheckScorer("SecondStageScorer",already_exist) ;
+
+ if(already_exist)
+ return ;
+
+ // Otherwise the scorer is initialised
+ G4VPrimitiveScorer* FirstStageScorer = new DSSDScorers::PS_Rectangle("FirstStageScorer",1,
+ TrapezoidBaseLarge,
+ TrapezoidHeight,
+ 128,128);
+ G4VPrimitiveScorer* SecondStageScorer = new DSSDScorers::PS_Rectangle("SecondStageScorer",1,
+ TrapezoidBaseLarge,
+ TrapezoidHeight,
+ 16,16);
+
+ G4VPrimitiveScorer* InteractionFirstStage = new InteractionScorers::PS_Interactions("InteractionFirstStage",ms_InterCoord,0);
+ G4VPrimitiveScorer* InteractionSecondStage = new InteractionScorers::PS_Interactions("InteractionSecondStage",ms_InterCoord,0);
+
+ // Register it to the multifunctionnal detector
+ m_FirstStageScorer->RegisterPrimitive(FirstStageScorer);
+ m_FirstStageScorer->RegisterPrimitive(InteractionFirstStage);
+ m_SecondStageScorer->RegisterPrimitive(SecondStageScorer);
+ m_SecondStageScorer->RegisterPrimitive(InteractionSecondStage);
+
+ G4SDManager::GetSDMpointer()->AddNewDetector(m_FirstStageScorer);
+ G4SDManager::GetSDMpointer()->AddNewDetector(m_SecondStageScorer);
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+////////////////////////////////////////////////////////////////////////////////
+// Construct Method to be pass to the DetectorFactory //
+////////////////////////////////////////////////////////////////////////////////
+NPS::VDetector* PISTA::Construct(){
+ return (NPS::VDetector*) new PISTA();
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+////////////////////////////////////////////////////////////////////////////////
+// Registering the construct method to the factory //
+////////////////////////////////////////////////////////////////////////////////
+extern"C" {
+ class proxy_nps_PISTA{
+ public:
+ proxy_nps_PISTA(){
+ NPS::DetectorFactory::getInstance()->AddToken("PISTA","PISTA");
+ NPS::DetectorFactory::getInstance()->AddDetector("PISTA",PISTA::Construct);
+ }
+ };
+
+ proxy_nps_PISTA p_nps_PISTA;
+}
diff --git a/NPSimulation/Detectors/PISTA/PISTA.hh b/NPSimulation/Detectors/PISTA/PISTA.hh
new file mode 100644
index 0000000000000000000000000000000000000000..a635ec357c58bb66172b1edf39ef53485c1742a9
--- /dev/null
+++ b/NPSimulation/Detectors/PISTA/PISTA.hh
@@ -0,0 +1,112 @@
+#ifndef PISTA_h
+#define PISTA_h 1
+/*****************************************************************************
+ * Copyright (C) 2009-2020 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: Pierre Morfouace contact address: pierre.morfouace2@cea.fr *
+ * *
+ * Creation Date : May 2020 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class describe PISTA simulation *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ *****************************************************************************/
+
+// C++ header
+#include <string>
+#include <vector>
+using namespace std;
+
+// G4 headers
+#include "G4ThreeVector.hh"
+#include "G4RotationMatrix.hh"
+#include "G4LogicalVolume.hh"
+#include "G4MultiFunctionalDetector.hh"
+
+// NPTool header
+#include "NPSVDetector.hh"
+#include "TPISTAData.h"
+#include "NPInputParser.h"
+
+class PISTA : public NPS::VDetector{
+ ////////////////////////////////////////////////////
+ /////// Default Constructor and Destructor /////////
+ ////////////////////////////////////////////////////
+ public:
+ PISTA() ;
+ virtual ~PISTA() ;
+
+ ////////////////////////////////////////////////////
+ /////// Specific Function of this Class ///////////
+ ////////////////////////////////////////////////////
+ public:
+ // Cartesian
+ void AddDetector(G4ThreeVector POS);
+ // Spherical
+ void AddDetector(double R,double Theta,double Phi);
+
+
+ G4LogicalVolume* BuildTrapezoidDetector();
+
+ private:
+ G4LogicalVolume* m_TrapezoidDetector;
+
+ ////////////////////////////////////////////////////
+ ////// Inherite from NPS::VDetector class /////////
+ ////////////////////////////////////////////////////
+ public:
+ // Read stream at Configfile to pick-up parameters of detector (Position,...)
+ // Called in DetecorConstruction::ReadDetextorConfiguration Method
+ void ReadConfiguration(NPL::InputParser) ;
+
+ // Construct detector and inialise sensitive part.
+ // Called After DetecorConstruction::AddDetector Method
+ void ConstructDetector(G4LogicalVolume* world) ;
+
+ // Add Detector branch to the EventTree.
+ // Called After DetecorConstruction::AddDetector Method
+ void InitializeRootOutput() ;
+
+ // Read sensitive part and fill the Root tree.
+ // Called at in the EventAction::EndOfEventAvtion
+ void ReadSensitive(const G4Event* event) ;
+
+ public: // Scorer
+ // Initialize all Scorer used by the MUST2Array
+ void InitializeScorers() ;
+
+ // Associated Scorer
+ G4MultiFunctionalDetector* m_FirstStageScorer ;
+ G4MultiFunctionalDetector* m_SecondStageScorer ;
+ ////////////////////////////////////////////////////
+ ///////////Event class to store Data////////////////
+ ////////////////////////////////////////////////////
+ private:
+ TPISTAData* m_Event ;
+
+ ////////////////////////////////////////////////////
+ ///////////////Private intern Data//////////////////
+ ////////////////////////////////////////////////////
+ private: // Geometry
+ // Detector Coordinate
+ vector<double> m_R;
+ vector<double> m_Theta;
+ vector<double> m_Phi;
+
+ // Visualisation Attribute
+ G4VisAttributes* m_VisTrap;
+
+ // Needed for dynamic loading of the library
+ public:
+ static NPS::VDetector* Construct();
+};
+#endif
diff --git a/NPSimulation/Process/BeamReaction.cc b/NPSimulation/Process/BeamReaction.cc
index dcfe2ce8292617ddf1ba611d454f14ef924d06a3..2a6dc19339a2265ebcf560f758c74e3859a00d17 100644
--- a/NPSimulation/Process/BeamReaction.cc
+++ b/NPSimulation/Process/BeamReaction.cc
@@ -320,7 +320,7 @@ void NPS::BeamReaction::DoIt(const G4FastTrack& fastTrack,
// Angles
// Shoot and Set a Random ThetaCM
m_Reaction.ShootRandomThetaCM();
- double phi = RandFlat::shoot() * 2. * pi;
+ double phi = G4RandFlat::shoot() * 2. * pi;
//////////////////////////////////////////////////
///// Momentum and angles from kinematics /////
@@ -453,8 +453,8 @@ void NPS::BeamReaction::DoIt(const G4FastTrack& fastTrack,
// Shoot and Set a Random ThetaCM
//m_QFS.ShootRandomThetaCM();
//m_QFS.ShootRandomPhiCM();
- double theta = RandFlat::shoot() * pi;
- double phi = RandFlat::shoot() * 2. * pi - pi; //rand in [-pi,pi]
+ double theta = G4RandFlat::shoot() * pi;
+ double phi = G4RandFlat::shoot() * 2. * pi - pi; //rand in [-pi,pi]
double momX = gRandom->Gaus(0.,m_QFS.GetMomentumSigma());
double momY = gRandom->Gaus(0.,m_QFS.GetMomentumSigma());
double momZ = gRandom->Gaus(0.,m_QFS.GetMomentumSigma());
diff --git a/NPSimulation/Simulation.cc b/NPSimulation/Simulation.cc
index ba0ce3916324317772413aaaf22dea1290a84b48..ed862a21df62f611f1149d2c3d1cc51f136af25e 100644
--- a/NPSimulation/Simulation.cc
+++ b/NPSimulation/Simulation.cc
@@ -67,6 +67,7 @@ int main(int argc, char** argv){
// initialize the state of the root and geant4 random generator
if(OptionManager->GetRandomSeed()>0){
+ std::cout << " Seeds for random generators set to: " << OptionManager->GetRandomSeed() << std::endl;
gRandom->SetSeed(OptionManager->GetRandomSeed());
CLHEP::HepRandom::setTheSeed(OptionManager->GetRandomSeed(),3);
}
diff --git a/Projects/Dali/PhysicsListOption.txt b/Projects/Dali/PhysicsListOption.txt
index a8991fd396fd70fc80b1cdcb14efc5fec098c56e..55051b2daa3ff7d02c80ce05ada9b18aae3fecb9 100644
--- a/Projects/Dali/PhysicsListOption.txt
+++ b/Projects/Dali/PhysicsListOption.txt
@@ -1,7 +1,7 @@
EmPhysicsList Option4
DefaultCutOff 10000
IonBinaryCascadePhysics 0
-DriftElectronPhysics 1
+DriftElectronPhysics 1
NPIonInelasticPhysics 0
EmExtraPhysics 0
HadronElasticPhysics 0
diff --git a/Projects/Minos/PhysicsListOption.txt b/Projects/Minos/PhysicsListOption.txt
index 5e948e6061566e703143f8c5a86e4ee5ff6d7f76..68f97e860e04654d462d3f6e1888f1ab694af8de 100644
--- a/Projects/Minos/PhysicsListOption.txt
+++ b/Projects/Minos/PhysicsListOption.txt
@@ -1,5 +1,5 @@
EmPhysicsList Option4
-DriftElectronPhysics 1
+DriftElectronPhysics 0
DefaultCutOff 1e9
IonBinaryCascadePhysics 0
NPIonInelasticPhysics 0
diff --git a/Projects/PISTA/Analysis.cxx b/Projects/PISTA/Analysis.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..ed8955b67a5196df8a1ff48574375837dddae3a2
--- /dev/null
+++ b/Projects/PISTA/Analysis.cxx
@@ -0,0 +1,164 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2016 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: XAUTHORX contact address: XMAILX *
+ * *
+ * Creation Date : XMONTHX XYEARX *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class describe PISTA analysis project *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ *****************************************************************************/
+
+#include<iostream>
+using namespace std;
+#include"Analysis.h"
+#include"NPAnalysisFactory.h"
+#include"NPDetectorManager.h"
+#include"NPOptionManager.h"
+#include"RootOutput.h"
+#include"RootInput.h"
+////////////////////////////////////////////////////////////////////////////////
+Analysis::Analysis(){
+}
+////////////////////////////////////////////////////////////////////////////////
+Analysis::~Analysis(){
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::Init(){
+ PISTA= (TPISTAPhysics*) m_DetectorManager->GetDetector("PISTA");
+ InitialConditions = new TInitialConditions();
+ ReactionConditions = new TReactionConditions();
+ InitOutputBranch();
+ InitInputBranch();
+ Rand = TRandom3();
+
+ TargetThickness = m_DetectorManager->GetTargetThickness();
+
+ Transfer = new NPL::Reaction("238U(12C,10Be)240Pu@1428");
+
+ // Energy loss table
+ Be10C = EnergyLoss("EnergyLossTable/Be10_C.G4table","G4Table",100);
+ U238C = EnergyLoss("EnergyLossTable/U238_C.G4table","G4Table",100);
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::TreatEvent(){
+ ReInitValue();
+
+ OriginalThetaLab = ReactionConditions->GetTheta(0);
+ OriginalElab = ReactionConditions->GetKineticEnergy(0);
+ OriginalBeamEnergy = ReactionConditions->GetBeamEnergy();
+
+ XTarget = InitialConditions->GetIncidentPositionX();
+ YTarget = InitialConditions->GetIncidentPositionY();
+ ZTarget = InitialConditions->GetIncidentPositionZ();
+
+ TVector3 BeamDirection = InitialConditions->GetBeamDirection();
+ TVector3 BeamPosition(XTarget,YTarget,ZTarget);
+
+ BeamEnergy = InitialConditions->GetIncidentInitialKineticEnergy();
+ BeamEnergy = U238C.Slow(BeamEnergy,TargetThickness*0.5,0);
+
+ Transfer->SetBeamEnergy(BeamEnergy);
+
+ for(unsigned int i = 0; i<PISTA->EventMultiplicity; i++){
+ if(PISTA->E.size()>0){
+ double Energy = PISTA->DE[i] + PISTA->E[i];
+
+ TVector3 HitDirection = PISTA->GetPositionOfInteraction(i);
+ //ThetaLab = HitDirection.Angle(BeamDirection);
+ ThetaLab = HitDirection.Angle(TVector3(0,0,1));
+
+ ThetaDetectorSurface = HitDirection.Angle(-PISTA->GetDetectorNormal(i));
+ ThetaNormalTarget = HitDirection.Angle(TVector3(0,0,1));
+
+ Elab = Be10C.EvaluateInitialEnergy(Energy,TargetThickness*0.5,ThetaNormalTarget);
+
+ OptimumEx = Transfer->ReconstructRelativistic(OriginalElab, OriginalThetaLab*deg);
+ Ex = Transfer->ReconstructRelativistic(Elab, ThetaLab);
+ ThetaCM = Transfer->EnergyLabToThetaCM(Elab, ThetaLab)/deg;
+ ThetaLab = ThetaLab/deg;
+ }
+ }
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::InitOutputBranch(){
+ RootOutput::getInstance()->GetTree()->Branch("OriginalBeamEnergy",&OriginalBeamEnergy,"OriginalBeamEnergy/D");
+ RootOutput::getInstance()->GetTree()->Branch("BeamEnergy",&BeamEnergy,"BeamEnergy/D");
+ RootOutput::getInstance()->GetTree()->Branch("XTarget",&XTarget,"XTarget/D");
+ RootOutput::getInstance()->GetTree()->Branch("YTarget",&YTarget,"YTarget/D");
+ RootOutput::getInstance()->GetTree()->Branch("ZTarget",&ZTarget,"ZTarget/D");
+ RootOutput::getInstance()->GetTree()->Branch("OptimumEx",&OptimumEx,"OptimumEx/D");
+ RootOutput::getInstance()->GetTree()->Branch("Ex",&Ex,"Ex/D");
+ RootOutput::getInstance()->GetTree()->Branch("Elab",&Elab,"Elab/D");
+ RootOutput::getInstance()->GetTree()->Branch("OriginalElab",&OriginalElab,"OriginalElab/D");
+ RootOutput::getInstance()->GetTree()->Branch("ThetaLab",&ThetaLab,"ThetaLab/D");
+ RootOutput::getInstance()->GetTree()->Branch("OriginalThetaLab",&OriginalThetaLab,"OriginalThetaLab/D");
+ RootOutput::getInstance()->GetTree()->Branch("ThetaCM",&ThetaCM,"ThetaCM/D");
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::InitInputBranch(){
+ RootInput::getInstance()->GetChain()->SetBranchStatus("InitialConditions",true);
+ RootInput::getInstance()->GetChain()->SetBranchStatus("fIC_*",true);
+ RootInput::getInstance()->GetChain()->SetBranchAddress("InitialConditions",&InitialConditions);
+ RootInput::getInstance()->GetChain()->SetBranchStatus("ReactionConditions",true);
+ RootInput::getInstance()->GetChain()->SetBranchStatus("fRC_*",true);
+ RootInput::getInstance()->GetChain()->SetBranchAddress("ReactionConditions",&ReactionConditions);
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::ReInitValue(){
+ OriginalBeamEnergy = -1000;
+ BeamEnergy = -1000;
+ OptimumEx = -1000;
+ Ex = -1000;
+ Elab = -1000;
+ OriginalElab = -1000;
+ OriginalThetaLab = -1000;
+ ThetaLab = -1000;
+ ThetaCM = -1000;
+ XTarget = -1000;
+ YTarget = -1000;
+ ZTarget = -1000;
+ OriginalThetaLab = -1000;
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::End(){
+}
+
+
+////////////////////////////////////////////////////////////////////////////////
+// Construct Method to be pass to the DetectorFactory //
+////////////////////////////////////////////////////////////////////////////////
+NPL::VAnalysis* Analysis::Construct(){
+ return (NPL::VAnalysis*) new Analysis();
+}
+
+////////////////////////////////////////////////////////////////////////////////
+// Registering the construct method to the factory //
+////////////////////////////////////////////////////////////////////////////////
+extern "C"{
+ class proxy{
+ public:
+ proxy(){
+ NPL::AnalysisFactory::getInstance()->SetConstructor(Analysis::Construct);
+ }
+ };
+
+ proxy p;
+}
+
diff --git a/Projects/PISTA/Analysis.h b/Projects/PISTA/Analysis.h
new file mode 100644
index 0000000000000000000000000000000000000000..3ec3cf61893310cc6e5a272963e05a2bee31b792
--- /dev/null
+++ b/Projects/PISTA/Analysis.h
@@ -0,0 +1,75 @@
+#ifndef Analysis_h
+#define Analysis_h
+/*****************************************************************************
+ * Copyright (C) 2009-2016 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: XAUTHORX contact address: XMAILX *
+ * *
+ * Creation Date : XMONTHX XYEARX *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class describe PISTA analysis project *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ *****************************************************************************/
+
+#include "NPVAnalysis.h"
+#include "TPISTAPhysics.h"
+#include "TInitialConditions.h"
+#include "TReactionConditions.h"
+#include "NPEnergyLoss.h"
+#include "NPReaction.h"
+#include "TRandom3.h"
+class Analysis: public NPL::VAnalysis{
+ public:
+ Analysis();
+ ~Analysis();
+
+ public:
+ void Init();
+ void TreatEvent();
+ void End();
+ void InitOutputBranch();
+ void InitInputBranch();
+ void ReInitValue();
+
+ static NPL::VAnalysis* Construct();
+
+ private:
+ double OriginalBeamEnergy;
+ double BeamEnergy;
+ double XTarget;
+ double YTarget;
+ double ZTarget;
+ double OriginalElab;
+ double Elab;
+ double ThetaLab;
+ double ThetaCM;
+ double OptimumEx;
+ double Ex;
+ double OriginalThetaLab;
+ NPL::Reaction* Transfer;
+
+ TRandom3 Rand;
+ double ThetaNormalTarget;
+ double ThetaDetectorSurface;
+ double TargetThickness;
+
+ NPL::EnergyLoss Be10C;
+ NPL::EnergyLoss U238C;
+
+ private:
+ TPISTAPhysics* PISTA;
+ TInitialConditions* InitialConditions;
+ TReactionConditions* ReactionConditions;
+
+};
+#endif
diff --git a/Projects/PISTA/CMakeLists.txt b/Projects/PISTA/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..22c74affdfc45019bdda2594f8439c52d4ab97ec
--- /dev/null
+++ b/Projects/PISTA/CMakeLists.txt
@@ -0,0 +1,5 @@
+cmake_minimum_required (VERSION 2.8)
+# 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/PISTA/PISTA.detector b/Projects/PISTA/PISTA.detector
new file mode 100644
index 0000000000000000000000000000000000000000..bc33042c2a7a66e9f543302264bf1f56b9298a04
--- /dev/null
+++ b/Projects/PISTA/PISTA.detector
@@ -0,0 +1,51 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Target
+ THICKNESS= 0.22 micrometer
+ RADIUS= 20 mm
+ MATERIAL= C
+ ANGLE= 0 deg
+ X= 0 mm
+ Y= 0 mm
+ Z= 0 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 315 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 270 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 225 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 180 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 135 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 90 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 45 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+PISTA
+ R= 100 mm
+ THETA= 60 deg
+ PHI= 0 deg
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
diff --git a/Projects/Scone/Analysis.cxx b/Projects/Scone/Analysis.cxx
index 6001665bb896524febfd3ae16bf16e189b094f84..c7210aba98a828545e57fd494f1091b8ddc22a07 100644
--- a/Projects/Scone/Analysis.cxx
+++ b/Projects/Scone/Analysis.cxx
@@ -84,9 +84,9 @@ void Analysis::TreatEvent(){
E_sum += Scone->Energy[i];
}
}
-
+ E_sum = E_sum - E_init;
//if(Time_max>50) m_DetectedNeutron++;
- if(Time_max>40 && E_sum>0.2) m_DetectedNeutron++;
+ if(Time_max>50 && E_sum>0) m_DetectedNeutron++;
}
////////////////////////////////////////////////////////////////////////////////
@@ -128,8 +128,9 @@ void Analysis::End(){
ofile.open("macro/eff_scone_natGd25um.txt");
//ofile.open("macro/eff_scone_menate.txt");
for(int i=0; i< vDetectedNeutron.size(); i++){
+ //cout << "* " << vE_init[i] << " / " << vDetectedNeutron[i]/vDetectedNeutron[0]*99.4 << endl;
cout << "* " << vE_init[i] << " / " << vDetectedNeutron[i]/1e5*100 << endl;
- //ofile << vE_init[i] << " " << vDetectedNeutron[i]/vDetectedNeutron[0]*99.3 << endl;
+ //ofile << vE_init[i] << " " << vDetectedNeutron[i]/vDetectedNeutron[0]*99.4 << endl;
ofile << vE_init[i] << " " << vDetectedNeutron[i]/1e5*100 << endl;
}
ofile.close();