diff --git a/NPLib/Detectors/ZDD/CMakeLists.txt b/NPLib/Detectors/ZDD/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..d902c986453e8c2f88dab6572856dd75558f69ae
--- /dev/null
+++ b/NPLib/Detectors/ZDD/CMakeLists.txt
@@ -0,0 +1,6 @@
+add_custom_command(OUTPUT TZDDPhysicsDict.cxx COMMAND ../../scripts/build_dict.sh TZDDPhysics.h TZDDPhysicsDict.cxx TZDDPhysics.rootmap libNPZDD.dylib DEPENDS TZDDPhysics.h)
+add_custom_command(OUTPUT TZDDDataDict.cxx COMMAND ../../scripts/build_dict.sh TZDDData.h TZDDDataDict.cxx TZDDData.rootmap libNPZDD.dylib DEPENDS TZDDData.h)
+add_library(NPZDD SHARED TZDDSpectra.cxx TZDDData.cxx TZDDPhysics.cxx TZDDDataDict.cxx TZDDPhysicsDict.cxx )
+target_link_libraries(NPZDD ${ROOT_LIBRARIES} NPCore)
+install(FILES TZDDData.h TZDDPhysics.h TZDDSpectra.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY})
+
diff --git a/NPLib/Detectors/ZDD/TZDDData.cxx b/NPLib/Detectors/ZDD/TZDDData.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..e1885d22c0e6ecdb26f2032e40a704571eb129de
--- /dev/null
+++ b/NPLib/Detectors/ZDD/TZDDData.cxx
@@ -0,0 +1,96 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold ZDD Raw data *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+#include "TZDDData.h"
+
+#include <iostream>
+#include <fstream>
+#include <sstream>
+#include <string>
+using namespace std;
+
+ClassImp(TZDDData)
+
+
+//////////////////////////////////////////////////////////////////////
+TZDDData::TZDDData() {
+}
+
+
+
+//////////////////////////////////////////////////////////////////////
+TZDDData::~TZDDData() {
+}
+
+
+
+//////////////////////////////////////////////////////////////////////
+void TZDDData::Clear() {
+ // Drift_Chambers
+ fZDD_Drift_DetectorNbr.clear();
+ fZDD_DriftTime.clear();
+ // Energy
+ fZDD_E_DetectorNbr.clear();
+ fZDD_Energy.clear();
+ // Time
+ fZDD_T_DetectorNbr.clear();
+ fZDD_Time.clear();
+
+ // IC Energy
+ fZDD_E_IC_Nbr.clear();
+ fZDD_IC_Energy.clear();
+ // IC Time
+ fZDD_T_IC_Nbr.clear();
+ fZDD_IC_Time.clear();
+
+ // Plastic Energy
+ fZDD_E_Plastic_Nbr.clear();
+ fZDD_Plastic_Energy.clear();
+ // Plastic Time
+ fZDD_T_Plastic_Nbr.clear();
+ fZDD_Plastic_Time.clear();
+}
+
+
+
+//////////////////////////////////////////////////////////////////////
+void TZDDData::Dump() const {
+ // This method is very useful for debuging and worth the dev.
+ cout << "XXXXXXXXXXXXXXXXXXXXXXXX New Event [TZDDData::Dump()] XXXXXXXXXXXXXXXXX" << endl;
+
+ // Energy
+ size_t mysize = fZDD_E_DetectorNbr.size();
+ cout << "ZDD_E_Mult: " << mysize << endl;
+
+ for (size_t i = 0 ; i < mysize ; i++){
+ cout << "DetNbr: " << fZDD_E_DetectorNbr[i]
+ << " Energy: " << fZDD_Energy[i];
+ }
+
+ // Time
+ mysize = fZDD_T_DetectorNbr.size();
+ cout << "ZDD_T_Mult: " << mysize << endl;
+
+ for (size_t i = 0 ; i < mysize ; i++){
+ cout << "DetNbr: " << fZDD_T_DetectorNbr[i]
+ << " Time: " << fZDD_Time[i];
+ }
+}
diff --git a/NPLib/Detectors/ZDD/TZDDData.h b/NPLib/Detectors/ZDD/TZDDData.h
new file mode 100644
index 0000000000000000000000000000000000000000..6ac5e4bfa6cc23ebcd775b49aa12c132ea483083
--- /dev/null
+++ b/NPLib/Detectors/ZDD/TZDDData.h
@@ -0,0 +1,219 @@
+#ifndef __ZDDDATA__
+#define __ZDDDATA__
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold ZDD Raw data *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+// STL
+#include <vector>
+#include <iostream>
+using namespace std;
+
+// ROOT
+#include "TObject.h"
+
+
+class TZDDData : public TObject {
+ //////////////////////////////////////////////////////////////
+ // data members are hold into vectors in order
+ // to allow multiplicity treatment
+ private:
+ double counter = 0;
+ // Energy IC
+ vector<UShort_t> fZDD_E_IC_Nbr;
+ vector<Double_t> fZDD_IC_Energy;
+
+ // Time IC
+ vector<UShort_t> fZDD_T_IC_Nbr;
+ vector<Double_t> fZDD_IC_Time;
+
+ // Energy Plastic
+ vector<UShort_t> fZDD_E_Plastic_Nbr;
+ vector<Double_t> fZDD_Plastic_Energy;
+
+ // Time Plastic
+ vector<UShort_t> fZDD_T_Plastic_Nbr;
+ vector<Double_t> fZDD_Plastic_Time;
+
+ // Energy
+ vector<UShort_t> fZDD_E_DetectorNbr;
+ vector<Double_t> fZDD_Energy;
+
+ // Time
+ vector<UShort_t> fZDD_T_DetectorNbr;
+ vector<Double_t> fZDD_Time;
+
+ // DriftTime in DC
+ vector<UShort_t> fZDD_Drift_DetectorNbr;
+ vector<Double_t> fZDD_DriftTime;
+
+ //////////////////////////////////////////////////////////////
+ // Constructor and destructor
+ public:
+ TZDDData();
+ ~TZDDData();
+
+
+ //////////////////////////////////////////////////////////////
+ // 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 ////////////////////////
+ // Energy
+ inline void SetEnergy(const UShort_t& DetNbr,const Double_t& Energy){
+ fZDD_E_DetectorNbr.push_back(DetNbr);
+ fZDD_Energy.push_back(Energy);
+ };//!
+
+ // Time
+ inline void SetTime(const UShort_t& DetNbr,const Double_t& Time) {
+ fZDD_T_DetectorNbr.push_back(DetNbr);
+ fZDD_Time.push_back(Time);
+ };//!
+
+ // IC Energy
+ inline void Set_IC_Energy(const UShort_t& IC_Nbr,const Double_t& IC_Energy){
+ fZDD_E_IC_Nbr.push_back(IC_Nbr);
+ fZDD_IC_Energy.push_back(IC_Energy);
+ };//!
+
+ // IC Time
+ inline void Set_IC_Time(const UShort_t& IC_Nbr,const Double_t& IC_Time) {
+ fZDD_T_IC_Nbr.push_back(IC_Nbr);
+ fZDD_IC_Time.push_back(IC_Time);
+ };//!
+
+ // Plastic Energy
+ inline void Set_Plastic_Energy(const UShort_t& Plastic_Nbr,const Double_t& Plastic_Energy){
+ fZDD_E_Plastic_Nbr.push_back(Plastic_Nbr);
+ fZDD_Plastic_Energy.push_back(Plastic_Energy);
+ };//!
+
+ // Plastic Time
+ inline void Set_Plastic_Time(const UShort_t& Plastic_Nbr,const Double_t& Plastic_Time) {
+ fZDD_T_Plastic_Nbr.push_back(Plastic_Nbr);
+ fZDD_Plastic_Time.push_back(Plastic_Time);
+ };//!
+
+ // Position DriftTime and X in DC
+ inline void Set_DC_Time(const UShort_t& DetNbr, const Double_t& DriftTime) {
+ fZDD_Drift_DetectorNbr.push_back(DetNbr);
+ fZDD_DriftTime.push_back(DriftTime);
+ }; //!
+
+ ////////////////////// GETTERS ////////////////////////
+ // Energy
+ inline UShort_t GetMultEnergy(std::string Detector) const
+ {if(Detector == "IC")
+ return fZDD_E_IC_Nbr.size();
+ else if(Detector == "Plastic")
+ return fZDD_E_Plastic_Nbr.size();
+ else{
+ std::cout << "Detector should be either IC or Plastic" << std::endl;
+ return -1;
+ }
+ }
+
+ inline UShort_t GetE_DetectorNbr(std::string Detector, const unsigned int &i) const
+ {if(Detector == "IC")
+ return fZDD_E_IC_Nbr[i];
+ else if(Detector == "Plastic")
+ return fZDD_E_Plastic_Nbr[i];
+ else{
+ std::cout << "Detector should be either IC or Plastic" << std::endl;
+ return -1;
+ }
+ }
+ inline Double_t Get_Energy(std::string Detector, const unsigned int &i) const
+ {if(Detector == "IC")
+ return fZDD_IC_Energy[i];
+ else if(Detector == "Plastic")
+ return fZDD_Plastic_Energy[i];
+ else{
+ std::cout << "Detector should be either IC or Plastic" << std::endl;
+ return -1;
+ }
+ }
+
+ // Time
+ inline UShort_t GetMultTime(std::string Detector) const
+ {if(Detector == "IC")
+ return fZDD_T_IC_Nbr.size();
+ else if(Detector == "Plastic")
+ return fZDD_T_Plastic_Nbr.size();
+ else if(Detector == "DC")
+ return fZDD_Drift_DetectorNbr.size();
+ else{
+ std::cout << "Detector should be either IC, DC or Plastic" << std::endl;
+ return -1;
+ }
+ }
+
+ inline UShort_t GetT_DetectorNbr(std::string Detector, const unsigned int &i) const
+ {if(Detector == "IC")
+ return fZDD_T_IC_Nbr[i];
+ else if(Detector == "Plastic")
+ return fZDD_T_Plastic_Nbr[i];
+ else if(Detector == "DC")
+ return fZDD_Drift_DetectorNbr[i];
+ else{
+ std::cout << "Detector should be either IC, DC or Plastic" << std::endl;
+ return -1;
+ }
+ }
+ inline Double_t Get_Time(std::string Detector, const unsigned int &i) const
+ {if(Detector == "IC")
+ return fZDD_IC_Time[i];
+ else if(Detector == "Plastic")
+ return fZDD_Plastic_Time[i];
+ else if(Detector == "DC")
+ return fZDD_DriftTime[i];
+ else{
+ std::cout << "Detector should be either IC, DC or Plastic" << std::endl;
+ return -1;
+ }
+ }
+ // Position
+ inline UShort_t GetMultDrift() const { return fZDD_Drift_DetectorNbr.size(); }
+ inline UShort_t GetDrift_DetectorNbr(const unsigned int& i) const {
+ return fZDD_Drift_DetectorNbr[i];
+ } //!
+ inline Double_t Get_DriftTime(const unsigned int& i) const {
+ return fZDD_DriftTime[i];
+ } //!
+
+
+ //////////////////////////////////////////////////////////////
+ // Required for ROOT dictionnary
+ ClassDef(TZDDData,1) // ZDDData structure
+};
+
+#endif
\ No newline at end of file
diff --git a/NPLib/Detectors/ZDD/TZDDPhysics.cxx b/NPLib/Detectors/ZDD/TZDDPhysics.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..ab4d6319ec1be05c79c52d9cc2c7e2faba5df00a
--- /dev/null
+++ b/NPLib/Detectors/ZDD/TZDDPhysics.cxx
@@ -0,0 +1,552 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold ZDD Treated data *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+#include "TZDDPhysics.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(TZDDPhysics)
+
+
+///////////////////////////////////////////////////////////////////////////
+TZDDPhysics::TZDDPhysics()
+ : m_EventData(new TZDDData),
+ m_PreTreatedData(new TZDDData),
+ m_EventPhysics(this),
+ m_Spectra(0),
+ m_E_RAW_Threshold(0), // adc channels
+ m_E_Threshold(0), // MeV
+ m_NumberOfDetectors(0)
+ {
+ ICcounter = 0;
+ ACcounter = 0;
+ GGcounter = 0;
+ Entry_Exit_counter = 0;
+ Plasticcounter = 0;
+}
+
+///////////////////////////////////////////////////////////////////////////
+/// A usefull method to bundle all operation to add a detector
+void TZDDPhysics::AddDetector(TVector3 , string ){
+ // 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 TZDDPhysics::AddDetector(double R, double Theta, double Phi, string shape){
+ // Compute the TVector3 corresponding
+ TVector3 Pos(R*sin(Theta)*cos(Phi),R*sin(Theta)*sin(Phi),R*cos(Theta));
+ // Call the cartesian method
+ AddDetector(Pos,shape);
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::BuildSimplePhysicalEvent() {
+ BuildPhysicalEvent();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::BuildPhysicalEvent() {
+ // apply thresholds and calibration
+ PreTreat();
+
+ // match energy and time together
+ Match_E_T("IC");
+ Match_E_T("Plastic");
+ Treat_DC();
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::Treat_DC(){
+ // Dont have anything to modify for the moment
+ unsigned int mysizeDC = m_PreTreatedData->GetMultDrift();
+ for(int i = 0; i < mysizeDC; i++){
+ DC_DetectorNumber.push_back(m_PreTreatedData->GetDrift_DetectorNbr(i));
+ DC_DriftTime.push_back(m_PreTreatedData->Get_DriftTime(i));
+ }
+}
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::Match_E_T(std::string Detector){
+ unsigned int mysizeE = m_PreTreatedData->GetMultEnergy(Detector);
+ unsigned int mysizeT = m_PreTreatedData->GetMultTime(Detector);
+ for (UShort_t e = 0; e < mysizeE ; e++) {
+ for (UShort_t t = 0; t < mysizeT ; t++) {
+ if (m_PreTreatedData->GetE_DetectorNbr(Detector, e) == m_PreTreatedData->GetT_DetectorNbr(Detector, t)) {
+ if(Detector == "IC"){
+ IC_DetectorNumber.push_back(m_PreTreatedData->GetE_DetectorNbr(Detector, e));
+ IC_Energy.push_back(m_PreTreatedData->Get_Energy(Detector, e));
+ IC_Time.push_back(m_PreTreatedData->Get_Time(Detector, t));
+ }
+ else if(Detector == "Plastic"){
+ Plastic_DetectorNumber.push_back(m_PreTreatedData->GetE_DetectorNbr(Detector, e));
+ Plastic_Energy.push_back(m_PreTreatedData->Get_Energy(Detector, e));
+ Plastic_Time.push_back(m_PreTreatedData->Get_Time(Detector, t));
+ }
+ else{
+ std::cout << "Detector should be IC or Plastic" << std::endl;
+ return;
+ }
+ }
+ }
+ }
+}
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::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();
+
+ // IC
+ PreTreatEnergy("IC", Cal);
+ PreTreatTime("IC", Cal);
+
+ // Plastic
+ PreTreatEnergy("Plastic", Cal);
+ PreTreatTime("Plastic", Cal);
+
+ //DC
+ PreTreatTime("DC", Cal);
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::PreTreatEnergy(std::string Detector, CalibrationManager* Cal){
+ unsigned int mysize = m_EventData->GetMultEnergy(Detector);
+ for (UShort_t i = 0; i < mysize ; ++i) {
+ if (m_EventData->Get_Energy(Detector, i) > m_E_RAW_Threshold) {
+ Double_t Energy = Cal->ApplyCalibration("ZDD/ENERGY"+NPL::itoa(m_EventData->GetE_DetectorNbr(Detector, i)),m_EventData->Get_Energy(Detector, i));
+ if (Energy > m_E_Threshold) {
+ if(Detector == "IC")
+ m_PreTreatedData->Set_IC_Energy(m_EventData->GetE_DetectorNbr(Detector, i), Energy);
+
+ else if(Detector == "Plastic")
+ m_PreTreatedData->Set_Plastic_Energy(m_EventData->GetE_DetectorNbr(Detector, i), Energy);
+ }
+ }
+ }
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::PreTreatTime(std::string Detector, CalibrationManager* Cal){
+ unsigned int mysize = m_EventData->GetMultTime(Detector);
+ for (UShort_t i = 0; i < mysize; ++i) {
+ Double_t Time= Cal->ApplyCalibration("ZDD/TIME"+NPL::itoa(m_EventData->GetT_DetectorNbr(Detector, i)),m_EventData->Get_Time(Detector,i));
+ if(Detector == "IC")
+ m_PreTreatedData->Set_IC_Time(m_EventData->GetT_DetectorNbr(Detector, i), Time);
+ else if(Detector == "Plastic")
+ m_PreTreatedData->Set_Plastic_Time(m_EventData->GetT_DetectorNbr(Detector, i), Time);
+ else if(Detector == "DC")
+ m_PreTreatedData->Set_DC_Time(m_EventData->GetT_DetectorNbr(Detector, i), Time);
+ }
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::ReadAnalysisConfig() {
+ bool ReadingStatus = false;
+
+ // path to file
+ string FileName = "./configs/ConfigZDD.dat";
+
+ // open analysis config file
+ ifstream AnalysisConfigFile;
+ AnalysisConfigFile.open(FileName.c_str());
+
+ if (!AnalysisConfigFile.is_open()) {
+ cout << " No ConfigZDD.dat found: Default parameter loaded for Analayis " << FileName << endl;
+ return;
+ }
+ cout << " Loading user parameter for Analysis from ConfigZDD.dat " << endl;
+
+ // Save it in a TAsciiFile
+ TAsciiFile* asciiConfig = RootOutput::getInstance()->GetAsciiFileAnalysisConfig();
+ asciiConfig->AppendLine("%%% ConfigZDD.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 = "ConfigZDD";
+ 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 TZDDPhysics::Clear() {
+ IC_DetectorNumber.clear();
+ IC_Energy.clear();
+ IC_Time.clear();
+ Plastic_DetectorNumber.clear();
+ Plastic_Energy.clear();
+ Plastic_Time.clear();
+ DC_DetectorNumber.clear();
+ DC_DriftTime.clear();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::ReadConfiguration(NPL::InputParser parser){
+
+ vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("ZDD");
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << "//// " << blocks.size() << " detectors found " << endl;
+
+ vector<string> TokenZDD = {"R", "Theta"};
+ vector<string> TokenDC = {"Z", "Gas","Thickness", "Pressure", "Temperature"};
+ vector<string> TokenAC = {"Z", "Thickness", "Material"};
+ vector<string> TokenEntryExit = {"Z", "Thickness", "Material"};
+ vector<string> TokenIC = {"Z", "Thickness", "Gas", "Pressure", "Temperature"};
+ vector<string> TokenPlastic = {"Material", "Width", "Length", "Thickness", "Pos"};
+
+
+
+ for(unsigned int i = 0 ; i < blocks.size() ; i++){
+ if (blocks[i]->HasTokenList(TokenZDD)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// ZDD " << i + 1 << endl;
+ double R = blocks[i]->GetDouble("R", "mm");
+ double Theta = blocks[i]->GetDouble("Theta", "deg");
+ Add_ZDD(R, Theta);
+ }
+ else if (blocks[i]->GetMainValue() == "DC"
+ && blocks[i]->HasTokenList(TokenDC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// DC" << i + 1 << endl;
+ double Z = blocks[i]->GetDouble("Z", "mm");
+ double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ string Gas = blocks[i]->GetString("Gas");
+ double Pressure = blocks[i]->GetDouble("Pressure", "bar");
+ double Temperature = blocks[i]->GetDouble("Temperature", "kelvin");
+ Add_Drift_Chamber(Z, Thickness, Gas, Pressure, Temperature);
+ }
+ else if (blocks[i]->GetMainValue() == "IC"
+ && blocks[i]->HasTokenList(TokenIC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// IC" << ICcounter+1 << endl;
+ double Z = blocks[i]->GetDouble("Z", "mm");
+ double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ string Gas = blocks[i]->GetString("Gas");
+ double Pressure = blocks[i]->GetDouble("Pressure", "bar");
+ double Temperature = blocks[i]->GetDouble("Temperature", "kelvin");
+ Add_Ionisation_Chamber(Z, Thickness, Gas, Pressure, Temperature);
+ ICcounter++;
+ }
+ else if (blocks[i]->GetMainValue() == "GasGap"
+ && blocks[i]->HasTokenList(TokenIC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// GasGap" << GGcounter+1 << endl;
+ double Z = blocks[i]->GetDouble("Z", "mm");
+ double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ string Gas = blocks[i]->GetString("Gas");
+ double Pressure = blocks[i]->GetDouble("Pressure", "bar");
+ double Temperature = blocks[i]->GetDouble("Temperature", "kelvin");
+ Add_Gas_Gap(Z, Thickness, Gas, Pressure, Temperature);
+ GGcounter++;
+ }
+ else if (blocks[i]->GetMainValue() == "AC"
+ && blocks[i]->HasTokenList(TokenAC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// AC" << ACcounter+1 << endl;
+ double Z = blocks[i]->GetDouble("Z", "mm");
+ double Thickness = blocks[i]->GetDouble("Thickness", "um");
+ string Material = blocks[i]->GetString("Material");
+ Add_AC(Z, Thickness, Material);
+ ACcounter++;
+ }
+ else if (blocks[i]->GetMainValue() == "EntryExit"
+ && blocks[i]->HasTokenList(TokenEntryExit)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// AC" << Entry_Exit_counter+1 << endl;
+ double Z = blocks[i]->GetDouble("Z", "mm");
+ double Thickness = blocks[i]->GetDouble("Thickness", "um");
+ string Material = blocks[i]->GetString("Material");
+ Add_Entry_Exit(Z, Thickness, Material);
+ Entry_Exit_counter++;
+ }
+ else if (blocks[i]->GetMainValue() == "Plastic"
+ && blocks[i]->HasTokenList(TokenPlastic)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// Plastic" << Plasticcounter + 1 << endl;
+ string Material = blocks[i]->GetString("Material");
+ double Width = blocks[i]->GetDouble("Width", "mm");
+ double Length = blocks[i]->GetDouble("Length", "mm");
+ double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ int x = (blocks[i]->GetTVector3("Pos", "mm")).X();
+ int y = (blocks[i]->GetTVector3("Pos", "mm")).Y();
+ int z = (blocks[i]->GetTVector3("Pos", "mm")).Z();
+ std::vector<int> Pos{x, y, z};
+ Add_Plastic(Material, Width, Length, Thickness, Pos);
+ Plasticcounter++;
+ }
+ else{
+ cout << "ERROR: check your input file formatting " << endl;
+ exit(1);
+ }
+ }
+ quicksort(&Entry_Exit_Z, &Entry_Exit_Thickness, &Entry_Exit_Material);
+ quicksort(&AC_Z, &AC_Thickness, &AC_Material);
+ quicksort(&m_Ionisation_Chamber_Z, &m_Ionisation_Chamber_Thickness, &m_Ionisation_Chamber_Gas, &m_Gas_Gap_Pressure, &m_Ionisation_Chamber_Temperature);
+ quicksort(&m_Gas_Gap_Z, &m_Gas_Gap_Thickness, &m_Gas_Gap_Gas, &m_Gas_Gap_Pressure, &m_Gas_Gap_Temperature);
+}
+
+void TZDDPhysics::quicksort(std::vector<double>* Z, std::vector<double>* Thickness, std::vector<string>* Material){
+ for(int i = 0; i < Z->size(); i++){
+ }
+ if(Z->size()>0){
+ double minZ;
+ int minindex;
+ for(int i = 0; i < Z->size()-1; i++){
+ minZ = (*Z)[i+1];
+ minindex = i+1;
+ for(int p = i; p < Z->size(); p++){
+ if((*Z)[p] < minZ){
+ minZ = (*Z)[p];
+ minindex = p;
+ }
+ }
+ double interZ = (*Z)[minindex];
+ (*Z)[minindex] = (*Z)[i];
+ (*Z)[i] = interZ;
+
+ double interThickness = (*Thickness)[minindex];
+ (*Thickness)[minindex] = (*Thickness)[i];
+ (*Thickness)[i] = interThickness;
+
+ string interMaterial = (*Material)[minindex];
+ (*Material)[minindex] = (*Material)[i];
+ (*Material)[i] = interMaterial;
+ }
+ }
+ for(int i = 0; i < Z->size(); i++){
+ }
+}
+
+
+void TZDDPhysics::quicksort(std::vector<double>* Z, std::vector<double>* Thickness, std::vector<string>* Gas, std::vector<double>* Pressure, std::vector<double>* Temperature){
+ for(int i = 0; i < Z->size(); i++){
+ }
+ if(Z->size()>0){
+ double minZ;
+ int minindex;
+ for(int i = 0; i < Z->size()-1; i++){
+ minZ = (*Z)[i+1];
+ minindex = i+1;
+ for(int p = i; p < Z->size(); p++){
+ if((*Z)[p] < minZ){
+ minZ = (*Z)[p];
+ minindex = p;
+ }
+ }
+ double interZ = (*Z)[minindex];
+ (*Z)[minindex] = (*Z)[i];
+ (*Z)[i] = interZ;
+
+ double interThickness = (*Thickness)[minindex];
+ (*Thickness)[minindex] = (*Thickness)[i];
+ (*Thickness)[i] = interThickness;
+
+ string interGas = (*Gas)[minindex];
+ (*Gas)[minindex] = (*Gas)[i];
+ (*Gas)[i] = interGas;
+
+ double interPressure = (*Pressure)[minindex];
+ (*Pressure)[minindex] = (*Pressure)[i];
+ (*Pressure)[i] = interPressure;
+
+ double interTemperature = (*Temperature)[minindex];
+ (*Temperature)[minindex] = (*Temperature)[i];
+ (*Temperature)[i] = interTemperature;
+ }
+ }
+ for(int i = 0; i < Z->size(); i++){
+ }
+}
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::InitSpectra() {
+ m_Spectra = new TZDDSpectra(m_NumberOfDetectors);
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::FillSpectra() {
+ m_Spectra -> FillRawSpectra(m_EventData);
+ m_Spectra -> FillPreTreatedSpectra(m_PreTreatedData);
+ m_Spectra -> FillPhysicsSpectra(m_EventPhysics);
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::CheckSpectra() {
+ m_Spectra->CheckSpectra();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::ClearSpectra() {
+ // To be done
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+map< string , TH1*> TZDDPhysics::GetSpectra() {
+ if(m_Spectra)
+ return m_Spectra->GetMapHisto();
+ else{
+ map< string , TH1*> empty;
+ return empty;
+ }
+}
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::WriteSpectra() {
+ m_Spectra->WriteSpectra();
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::AddParameterToCalibrationManager() {
+ CalibrationManager* Cal = CalibrationManager::getInstance();
+ for (int i = 0; i < m_NumberOfDetectors; ++i) {
+ Cal->AddParameter("ZDD", "D"+ NPL::itoa(i+1)+"_ENERGY","ZDD_D"+ NPL::itoa(i+1)+"_ENERGY");
+ Cal->AddParameter("ZDD", "D"+ NPL::itoa(i+1)+"_TIME","ZDD_D"+ NPL::itoa(i+1)+"_TIME");
+ }
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::InitializeRootInputRaw() {
+ TChain* inputChain = RootInput::getInstance()->GetChain();
+ inputChain->SetBranchStatus("ZDD", true );
+ inputChain->SetBranchAddress("ZDD", &m_EventData );
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::InitializeRootInputPhysics() {
+ TChain* inputChain = RootInput::getInstance()->GetChain();
+ inputChain->SetBranchAddress("ZDD", &m_EventPhysics);
+}
+
+
+
+///////////////////////////////////////////////////////////////////////////
+void TZDDPhysics::InitializeRootOutput() {
+ TTree* outputTree = RootOutput::getInstance()->GetTree();
+ outputTree->Branch("ZDD", "TZDDPhysics", &m_EventPhysics);
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+// Construct Method to be pass to the DetectorFactory //
+////////////////////////////////////////////////////////////////////////////////
+NPL::VDetector* TZDDPhysics::Construct() {
+ return (NPL::VDetector*) new TZDDPhysics();
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+// Registering the construct method to the factory //
+////////////////////////////////////////////////////////////////////////////////
+extern "C"{
+class proxy_ZDD{
+ public:
+ proxy_ZDD(){
+ NPL::DetectorFactory::getInstance()->AddToken("ZDD","ZDD");
+ NPL::DetectorFactory::getInstance()->AddDetector("ZDD",TZDDPhysics::Construct);
+ }
+};
+
+proxy_ZDD p_ZDD;
+}
+
diff --git a/NPLib/Detectors/ZDD/TZDDPhysics.h b/NPLib/Detectors/ZDD/TZDDPhysics.h
new file mode 100644
index 0000000000000000000000000000000000000000..c9c32a74b291d9fc8e625cc260c0ab3f2aa4accf
--- /dev/null
+++ b/NPLib/Detectors/ZDD/TZDDPhysics.h
@@ -0,0 +1,335 @@
+#ifndef TZDDPHYSICS_H
+#define TZDDPHYSICS_H
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold ZDD 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 "TZDDData.h"
+#include "TZDDSpectra.h"
+#include "NPCalibrationManager.h"
+#include "NPVDetector.h"
+#include "NPInputParser.h"
+// forward declaration
+class TZDDSpectra;
+
+
+
+class TZDDPhysics : public TObject, public NPL::VDetector {
+ //////////////////////////////////////////////////////////////
+ // constructor and destructor
+ public:
+ TZDDPhysics();
+ ~TZDDPhysics() {};
+
+
+ //////////////////////////////////////////////////////////////
+ // Inherited from TObject and overriden to avoid warnings
+ public:
+ void Clear();
+ void Clear(const Option_t*) {};
+
+
+ //////////////////////////////////////////////////////////////
+ // data obtained after BuildPhysicalEvent() and stored in
+ // output ROOT file
+ public:
+ vector<int> IC_DetectorNumber;
+ vector<double> IC_Energy;
+ vector<double> IC_Time;
+
+ vector<int> Plastic_DetectorNumber;
+ vector<double> Plastic_Energy;
+ vector<double> Plastic_Time;
+
+ vector<int> DC_DetectorNumber;
+ vector<double> DC_DriftTime;
+
+ /// A usefull method to bundle all operation to add a detector
+ void AddDetector(TVector3 POS, string shape);
+ void AddDetector(double R, double Theta, double Phi, string shape);
+
+ //////////////////////////////////////////////////////////////
+ // methods inherited from the VDetector ABC class
+ public:
+ // read stream from ConfigFile to pick-up detector parameters
+ void ReadConfiguration(NPL::InputParser);
+
+ // 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 TZDDSpectra class
+ // instantiate the TZDDSpectra 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();
+
+ void Add_ZDD(double R, double Theta) {
+ /* cout << "Balise 2" << endl;
+ m_R = R;
+ m_Theta = Theta;
+ */
+ m_NumberOfDetectors++;
+ }
+
+ void Add_Ionisation_Chamber(double Z, double Thickness, string Gas, double Pressure,
+ double Temperature) {
+ m_Ionisation_Chamber_Z.push_back(Z);
+ m_Ionisation_Chamber_Thickness.push_back(Thickness);
+ m_Ionisation_Chamber_Gas.push_back(Gas);
+ m_Ionisation_Chamber_Pressure.push_back(Pressure);
+ m_Ionisation_Chamber_Temperature.push_back(Temperature);
+
+ m_NumberOfDetectors++;
+ }
+
+ void Add_Plastic(string Material, double Width, double Length,
+ double Thickness, std::vector<int> Pos) {
+/* m_Plastic_Material.push_back(Material);
+ m_Plastic_Width.push_back(Width);
+ m_Plastic_Length.push_back(Length);
+ m_Plastic_Thickness.push_back(Thickness);
+ m_Plastic_Position.push_back(Pos);
+ */
+ m_NumberOfDetectors++;
+ }
+ void Add_Gas_Gap(double Z,double Thickness,string Gas,double Pressure,double Temperature){
+ m_Gas_Gap_Z.push_back(Z);
+ m_Gas_Gap_Thickness.push_back(Thickness);
+ m_Gas_Gap_Gas.push_back(Gas);
+ m_Gas_Gap_Pressure.push_back(Pressure);
+ m_Gas_Gap_Temperature.push_back(Temperature);
+ m_NumberOfDetectors++;
+ };
+
+ void Add_AC(double Z, double Thickness, string Material){
+ AC_Material.push_back(Material);
+ AC_Thickness.push_back(Thickness);
+ AC_Z.push_back(Z);
+ m_NumberOfDetectors++;
+ };
+
+ void Add_Entry_Exit(double Z, double Thickness, string Material){
+ Entry_Exit_Material.push_back(Material);
+ Entry_Exit_Thickness.push_back(Thickness);
+ Entry_Exit_Z.push_back(Z);
+ m_NumberOfDetectors++;
+ };
+
+ void Add_Drift_Chamber(double Z, double Thickness, string Gas, double Pressure,
+ double Temperature) {
+ m_Drift_Chamber_Z.push_back(Z);
+ m_Drift_Chamber_Thickness.push_back(Thickness);
+ m_Drift_Chamber_Gas.push_back(Gas);
+ m_Drift_Chamber_Pressure.push_back(Pressure);
+ m_Drift_Chamber_Temperature.push_back(Temperature);
+
+ m_NumberOfDetectors++;
+ }
+
+
+ //////////////////////////////////////////////////////////////
+ // specific methods to ZDD array
+ public:
+ // remove bad channels, calibrate the data and apply thresholds
+ void PreTreat();
+
+
+ void Treat_DC();
+ // Matching Time and E for IC and Plastic
+ void Match_E_T(std::string Detector);
+
+ // PreTreating Energy for IC and Plastic
+ void PreTreatEnergy(std::string Detector, CalibrationManager* Cal);
+
+ // Same for time
+ void PreTreatTime(std::string Detector, CalibrationManager* Cal);
+
+ // 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 TZDDData object to TZDDPhysics.
+ // needed for online analysis for example
+ void SetRawDataPointer(TZDDData* rawDataPointer) {m_EventData = rawDataPointer;}
+
+ void quicksort(std::vector<double>* Z, std::vector<double>* Thickness, std::vector<string>* Material);
+ void quicksort(std::vector<double>* Z, std::vector<double>* Thickness, std::vector<string>* Gas, std::vector<double>* Pressure, std::vector<double>* Temperature);
+
+ // objects are not written in the TTree
+ private:
+ TZDDData* m_EventData; //!
+ TZDDData* m_PreTreatedData; //!
+ TZDDPhysics* m_EventPhysics; //!
+
+ // getters for raw and pre-treated data object
+ public:
+ TZDDData* GetRawData() const {return m_EventData;}
+ TZDDData* GetPreTreatedData() const {return m_PreTreatedData;}
+
+ // parameters used in the analysis
+ private:
+ // thresholds
+ double m_E_RAW_Threshold; //!
+ double m_E_Threshold; //!
+
+ // Nb of detectors counter
+ int ICcounter;
+ int ACcounter;
+ int GGcounter;
+ int Entry_Exit_counter;
+ int Plasticcounter;
+
+ public:
+ int GetICcounter(){return ICcounter;};
+
+ // number of detectors
+ private:
+ int m_NumberOfDetectors; //!
+ std::vector<string> Entry_Exit_Material; //!
+ std::vector<double> Entry_Exit_Thickness; //!
+ std::vector<double> Entry_Exit_Z; //!
+
+ std::vector<string> AC_Material; //!
+ std::vector<double> AC_Thickness; //!
+ std::vector<double> AC_Z; //!
+
+ std::vector<double> m_Ionisation_Chamber_Z; //!
+ std::vector<double> m_Ionisation_Chamber_Thickness; //!
+ std::vector<string> m_Ionisation_Chamber_Gas; //!
+ std::vector<double> m_Ionisation_Chamber_Pressure; //!
+ std::vector<double> m_Ionisation_Chamber_Temperature; //!
+
+ std::vector<double> m_Drift_Chamber_Z; //!
+ std::vector<double> m_Drift_Chamber_Thickness; //!
+ std::vector<string> m_Drift_Chamber_Gas; //!
+ std::vector<double> m_Drift_Chamber_Pressure; //!
+ std::vector<double> m_Drift_Chamber_Temperature; //!
+
+ std::vector<double> m_Gas_Gap_Z; //!
+ std::vector<double> m_Gas_Gap_Thickness; //!
+ std::vector<string> m_Gas_Gap_Gas; //!
+ std::vector<double> m_Gas_Gap_Pressure; //!
+ std::vector<double> m_Gas_Gap_Temperature; //!
+
+/* std::vector<double> CATS_X_1;
+ std::vector<double> CATS_X_2;
+ std::vector<double> CATS_Y_1;
+ std::vector<double> CATS_Y_2;
+
+ std::vector<double> CATS_TS_1;
+ std::vector<double> CATS_TS_2;*/
+
+ public:
+ std::vector<string> Get_Entry_Exit_Material(){return Entry_Exit_Material;};
+ std::vector<double> Get_Entry_Exit_Thickness(){return Entry_Exit_Thickness;};
+ std::vector<double> Get_Entry_Exit_Z(){return Entry_Exit_Z;};
+
+ std::vector<string> Get_AC_Material(){return AC_Material;};
+ std::vector<double> Get_AC_Thickness(){return AC_Thickness;};
+ std::vector<double> Get_AC_Z(){return AC_Z;};
+
+ std::vector<double>Get_m_Ionisation_Chamber_Z(){return m_Ionisation_Chamber_Z;};
+ std::vector<double>Get_m_Ionisation_Chamber_Thickness(){return m_Ionisation_Chamber_Thickness;};
+ std::vector<string>Get_m_Ionisation_Chamber_Gas(){return m_Ionisation_Chamber_Gas;};
+ std::vector<double>Get_m_Ionisation_Chamber_Pressure(){return m_Ionisation_Chamber_Pressure;};
+ std::vector<double>Get_m_Ionisation_Chamber_Temperature(){return m_Ionisation_Chamber_Temperature;};
+
+ std::vector<double>Get_m_Gas_Gap_Z(){return m_Gas_Gap_Z;};
+ std::vector<double>Get_m_Gas_Gap_Thickness(){return m_Gas_Gap_Thickness;};
+ std::vector<string>Get_m_Gas_Gap_Gas(){return m_Gas_Gap_Gas;};
+ std::vector<double>Get_m_Gas_Gap_Pressure(){return m_Gas_Gap_Pressure;};
+ std::vector<double>Get_m_Gas_Gap_Temperature(){return m_Gas_Gap_Temperature;};
+
+ std::vector<double>Get_m_Drift_Chamber_Z(){return m_Drift_Chamber_Z;};
+ std::vector<double>Get_m_Drift_Chamber_Thickness(){return m_Drift_Chamber_Thickness;};
+ std::vector<string>Get_m_Drift_Chamber_Gas(){return m_Drift_Chamber_Gas;};
+ std::vector<double>Get_m_Drift_Chamber_Pressure(){return m_Drift_Chamber_Pressure;};
+ std::vector<double>Get_m_Drift_Chamber_Temperature(){return m_Drift_Chamber_Temperature;};
+
+
+ // spectra class
+ private:
+ TZDDSpectra* m_Spectra; // !
+
+ // spectra getter
+ public:
+ map<string, TH1*> GetSpectra();
+
+ // Static constructor to be passed to the Detector Factory
+ public:
+ static NPL::VDetector* Construct();
+
+ ClassDef(TZDDPhysics,1) // ZDDPhysics structure
+};
+#endif
diff --git a/NPLib/Detectors/ZDD/TZDDSpectra.cxx b/NPLib/Detectors/ZDD/TZDDSpectra.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..37e9f081a9440eaf119b0eccd25cd6bcf7219fcb
--- /dev/null
+++ b/NPLib/Detectors/ZDD/TZDDSpectra.cxx
@@ -0,0 +1,174 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold ZDD Spectra *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+// class header
+#include "TZDDSpectra.h"
+
+// STL
+#include <iostream>
+#include <string>
+using namespace std;
+
+// NPTool header
+#include "NPOptionManager.h"
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+TZDDSpectra::TZDDSpectra()
+ : fNumberOfDetectors(0) {
+ SetName("ZDD");
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+TZDDSpectra::TZDDSpectra(unsigned int NumberOfDetectors) {
+ if(NPOptionManager::getInstance()->GetVerboseLevel()>0)
+ cout << "************************************************" << endl
+ << "TZDDSpectra : Initalizing control spectra for "
+ << NumberOfDetectors << " Detectors" << endl
+ << "************************************************" << endl ;
+ SetName("ZDD");
+ fNumberOfDetectors = NumberOfDetectors;
+
+ InitRawSpectra();
+ InitPreTreatedSpectra();
+ InitPhysicsSpectra();
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+TZDDSpectra::~TZDDSpectra() {
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TZDDSpectra::InitRawSpectra() {
+ static string name;
+ for (unsigned int i = 0; i < fNumberOfDetectors; i++) { // loop on number of detectors
+ // Energy
+ name = "ZDD"+NPL::itoa(i+1)+"_ENERGY_RAW";
+ AddHisto1D(name, name, 4096, 0, 16384, "ZDD/RAW");
+ // Time
+ name = "ZDD"+NPL::itoa(i+1)+"_TIME_RAW";
+ AddHisto1D(name, name, 4096, 0, 16384, "ZDD/RAW");
+ } // end loop on number of detectors
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TZDDSpectra::InitPreTreatedSpectra() {
+ static string name;
+ for (unsigned int i = 0; i < fNumberOfDetectors; i++) { // loop on number of detectors
+ // Energy
+ name = "ZDD"+NPL::itoa(i+1)+"_ENERGY_CAL";
+ AddHisto1D(name, name, 500, 0, 25, "ZDD/CAL");
+ // Time
+ name = "ZDD"+NPL::itoa(i+1)+"_TIME_CAL";
+ AddHisto1D(name, name, 500, 0, 25, "ZDD/CAL");
+
+
+ } // end loop on number of detectors
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TZDDSpectra::InitPhysicsSpectra() {
+ static string name;
+ // Kinematic Plot
+ name = "ZDD_ENERGY_TIME";
+ AddHisto2D(name, name, 500, 0, 500, 500, 0, 50, "ZDD/PHY");
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TZDDSpectra::FillRawSpectra(TZDDData* RawData) {
+ static string name;
+ static string family;
+
+ // Energy
+ unsigned int sizeE = RawData->GetMultEnergy("Plastic");
+ for (unsigned int i = 0; i < sizeE; i++) {
+ name = "ZDD"+NPL::itoa(RawData->GetE_DetectorNbr("Plastic",i))+"_ENERGY_RAW";
+ family = "ZDD/RAW";
+
+ FillSpectra(family,name,RawData->Get_Energy("Plastic", i));
+ }
+
+ // Time
+ unsigned int sizeT = RawData->GetMultTime("Plastic");
+ for (unsigned int i = 0; i < sizeT; i++) {
+ name = "ZDD"+NPL::itoa(RawData->GetT_DetectorNbr("Plastic", i))+"_TIME_RAW";
+ family = "ZDD/RAW";
+
+ FillSpectra(family,name,RawData->Get_Time("Plastic", i));
+ }
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TZDDSpectra::FillPreTreatedSpectra(TZDDData* PreTreatedData) {
+ static string name;
+ static string family;
+
+ // Energy
+ unsigned int sizeE = PreTreatedData->GetMultEnergy("Plastic");
+ for (unsigned int i = 0; i < sizeE; i++) {
+ name = "ZDD"+NPL::itoa(PreTreatedData->GetE_DetectorNbr("Plastic", i))+"_ENERGY_CAL";
+ family = "ZDD/CAL";
+
+ FillSpectra(family,name,PreTreatedData->Get_Energy("Plastic", i));
+ }
+
+ // Time
+ unsigned int sizeT = PreTreatedData->GetMultTime("Plastic");
+ for (unsigned int i = 0; i < sizeT; i++) {
+ name = "ZDD"+NPL::itoa(PreTreatedData->GetT_DetectorNbr("Plastic", i))+"_TIME_CAL";
+ family = "ZDD/CAL";
+
+ FillSpectra(family,name,PreTreatedData->Get_Time("Plastic", i));
+ }
+}
+
+
+
+////////////////////////////////////////////////////////////////////////////////
+void TZDDSpectra::FillPhysicsSpectra(TZDDPhysics* Physics) {
+ static string name;
+ static string family;
+ family= "ZDD/PHY";
+
+ // Energy vs time
+ unsigned int sizeE = Physics->Plastic_Energy.size();
+ for(unsigned int i = 0 ; i < sizeE ; i++){
+ name = "ZDD_ENERGY_TIME";
+ FillSpectra(family,name,Physics->Plastic_Energy[i],Physics->Plastic_Time[i]);
+ }
+}
+
diff --git a/NPLib/Detectors/ZDD/TZDDSpectra.h b/NPLib/Detectors/ZDD/TZDDSpectra.h
new file mode 100644
index 0000000000000000000000000000000000000000..925ebe8b540dc72f1cf142eb866fad6340b47cd7
--- /dev/null
+++ b/NPLib/Detectors/ZDD/TZDDSpectra.h
@@ -0,0 +1,62 @@
+#ifndef TZDDSPECTRA_H
+#define TZDDSPECTRA_H
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class hold ZDD Spectra *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ * *
+ *****************************************************************************/
+
+// NPLib headers
+#include "NPVSpectra.h"
+#include "TZDDData.h"
+#include "TZDDPhysics.h"
+
+// Forward Declaration
+class TZDDPhysics;
+
+
+class TZDDSpectra : public VSpectra {
+ //////////////////////////////////////////////////////////////
+ // constructor and destructor
+ public:
+ TZDDSpectra();
+ TZDDSpectra(unsigned int NumberOfDetectors);
+ ~TZDDSpectra();
+
+ //////////////////////////////////////////////////////////////
+ // Initialization methods
+ private:
+ void InitRawSpectra();
+ void InitPreTreatedSpectra();
+ void InitPhysicsSpectra();
+
+ //////////////////////////////////////////////////////////////
+ // Filling methods
+ public:
+ void FillRawSpectra(TZDDData*);
+ void FillPreTreatedSpectra(TZDDData*);
+ void FillPhysicsSpectra(TZDDPhysics*);
+
+ //////////////////////////////////////////////////////////////
+ // Detector parameters
+ private:
+ unsigned int fNumberOfDetectors;
+};
+
+#endif
diff --git a/NPSimulation/Detectors/ZDD/CMakeLists.txt b/NPSimulation/Detectors/ZDD/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..df95c2f9faa3e7c41e31fccf75d9bdd20daafde0
--- /dev/null
+++ b/NPSimulation/Detectors/ZDD/CMakeLists.txt
@@ -0,0 +1,2 @@
+add_library(NPSZDD SHARED ZDD.cc)
+target_link_libraries(NPSZDD NPSCore ${ROOT_LIBRARIES} ${Geant4_LIBRARIES} ${NPLib_LIBRARIES} -lNPZDD)
diff --git a/NPSimulation/Detectors/ZDD/ZDD.cc b/NPSimulation/Detectors/ZDD/ZDD.cc
new file mode 100644
index 0000000000000000000000000000000000000000..b928c6ef40c0b12a6ab148ba3ee41529fecbde10
--- /dev/null
+++ b/NPSimulation/Detectors/ZDD/ZDD.cc
@@ -0,0 +1,757 @@
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class describe ZDD simulation *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ *****************************************************************************/
+
+// C++ headers
+#include <sstream>
+#include <cmath>
+#include <limits>
+//G4 Geometry object
+#include "G4Tubs.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 "ZDD.hh"
+#include "DriftChamberScorers.hh"
+#include "CalorimeterScorers.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 ZDD_NS{
+ // Energy and time Resolution
+ /* const double EnergyThreshold = 0.1*MeV;
+ const double ResoTime = 4.5*ns ;
+ const double ResoEnergy = 1.0*MeV ;
+ const double Radius = 50*mm ;
+ const double Width = 100*mm ;
+ const double Thickness = 300*mm ;
+ const string Material = "BC400";
+ */
+ // IC E and T
+ const double EnergyThreshold_IC = 0.1*MeV;
+ const double ResoTime_IC = 4.5*ns ;
+ const double ResoEnergy_IC = 4; //Resolution en pourcentage
+
+ // Plastic E and T
+ const double EnergyThreshold_Plastic = 0.1*keV;
+ const double ResoTime_Plastic = 4.5*ns ;
+ const double ResoEnergy_Plastic = 5; // Resolution en pourcentage
+
+ // Drift Time and position of drift chambers
+ const G4double ResoDriftTime = 0.0001 * ns;
+
+ // Drift features
+ const G4double DriftSpeed = 1 * cm /microsecond; //microsecond; RQ très importante: le temps enregistré par les DC est un temps interne, ie il prend uniquement en compte le temps de dérive dans les drift
+ // En sortant des scorers, ce temps est en us
+ const G4ThreeVector DriftDir1 = G4ThreeVector(0, 1, 0);
+ const G4ThreeVector DriftDir2 = G4ThreeVector(1, 0, 0);
+
+
+ // Drift Chambers
+ const G4double Drift_Chamber_Width = 400 * mm;
+ const G4double Drift_Chamber_Length = 400 * mm;
+ const G4double Drift_Chamber_Thickness = 120 * mm;
+
+ // ZDD Volume
+ const G4double Phi = 0 * deg;
+ const G4double ZDD_Volume_Width = 4500 * mm;
+ const G4double ZDD_Volume_Length = 1000 * mm;
+ const G4double ZDD_Volume_Thickness = 5000 * mm;
+
+ // Anodes and Cathodes in IC (Mylar)
+ const G4double AC_Width = 400 * mm;
+ const G4double AC_Length = 400 * mm;
+// const G4double AC_Thickness = 2.31 * um;
+ //const G4double AC_Angle;
+
+ // IC Entrance and Exit (Kapton)
+ const G4double Kapton_Foil_Width = 400 * mm;
+ const G4double Kapton_Foil_Length = 400 * mm;
+// const G4doublKaptonar_Foil_Thickness = 2.31 * um;
+ const G4double Kapton_Foil_Angle = 30 * deg;
+
+ // Ionisation Chambers
+ const G4double Ionisation_Chamber_Width = 400 * mm;
+ const G4double Ionisation_Chamber_Length = 400 * mm;
+ //const G4double Ionisation_Chamber_Thickness = 40 * mm;
+ const G4double Ionisation_Chamber_Angle = 30 * deg;
+
+ // Gas Gap
+ const G4double Gas_Gap_Width = 400 * mm;
+ const G4double Gas_Gap_Length = 400 * mm;
+ //const G4double Gas_Gap_Thickness = 40 * mm;
+ const G4double Gas_Gap_Angle = 30 * deg;
+}
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+// ZDD Specific Method
+ZDD::ZDD(){
+ m_Event = new TZDDData() ;
+ m_ZDDScorer = 0;
+
+ ClearGeometry();
+ m_IC_Scorer = 0;
+ m_Plastic_Scorer = 0;
+ m_Drift_Chamber_Scorer_1 = 0;
+ m_Drift_Chamber_Scorer_2 = 0;
+ m_Drift_Chamber_1 = 0;
+ m_Drift_Chamber_2 = 0;
+
+ ICcounter = 0;
+ ACcounter = 0;
+ GGcounter = 0;
+ Entry_Exit_counter = 0;
+ Plasticcounter = 0;
+
+
+ // RGB Color + Transparency
+ m_VisSquare = new G4VisAttributes(G4Colour(0, 1, 0, 0.5));
+ m_VisCylinder = new G4VisAttributes(G4Colour(0, 0, 1, 0.5));
+ m_Vis_ZDD = new G4VisAttributes(G4Colour(1.0, 0.0, 0.0, 0.1));
+ m_Vis_Drift_Chamber_Gas = new G4VisAttributes(G4Colour(0, 1, 0, 0.5));
+ m_Vis_Ionisation_Chamber_Gas = new G4VisAttributes(G4Colour(0.0, 0.0, 1.0, 0.3));
+ m_Vis_AC = new G4VisAttributes(G4Colour(0.0, 0.5, 0.5, 0.5));
+ m_Vis_Plastic = new G4VisAttributes(G4Colour(1.0, 0.0, 0.0, 0.8));
+
+}
+
+ZDD::~ZDD(){
+}
+
+using namespace ZDD_NS;
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+/*void ZDD::AddDetector(G4ThreeVector POS, string Shape){
+ // 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());
+ m_Shape.push_back(Shape);
+}
+
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+void ZDD::AddDetector(double R, double Theta, double Phi, string Shape){
+ m_R.push_back(R);
+ m_Theta.push_back(Theta);
+ m_Phi.push_back(Phi);
+ m_Shape.push_back(Shape);
+}
+*/
+void ZDD::Add_ZDD(G4double R, double Theta) {
+ m_R = R;
+ m_Theta = Theta;
+}
+void ZDD::Add_Drift_Chamber(G4double Z, double Thickness, string Gas, double Pressure,double Temperature) {
+ m_Drift_Chamber_Z.push_back(Z);
+ m_Drift_Chamber_Gas.push_back(Gas);
+ m_Drift_Chamber_Pressure.push_back(Pressure);
+ m_Drift_Chamber_Temperature.push_back(Temperature);
+ m_Drift_Chamber_Thickness.push_back(Thickness);
+}
+
+void ZDD::Add_Ionisation_Chamber(G4double Z, double Thickness, string Gas, double Pressure,
+ double Temperature) {
+ m_Ionisation_Chamber_Z.push_back(Z);
+ m_Ionisation_Chamber_Thickness.push_back(Thickness);
+ m_Ionisation_Chamber_Gas.push_back(Gas);
+ m_Ionisation_Chamber_Pressure.push_back(Pressure);
+ m_Ionisation_Chamber_Temperature.push_back(Temperature);
+}
+
+void ZDD::Add_Gas_Gap(G4double Z, double Thickness, string Gas, double Pressure,
+ double Temperature) {
+ m_Gas_Gap_Z.push_back(Z);
+ m_Gas_Gap_Thickness.push_back(Thickness);
+ m_Gas_Gap_Gas.push_back(Gas);
+ m_Gas_Gap_Pressure.push_back(Pressure);
+ m_Gas_Gap_Temperature.push_back(Temperature);
+}
+
+void ZDD::Add_AC(G4double Z, double Thickness, string Material, G4double Angle) {
+ m_AC_Z.push_back(Z);
+ m_AC_Thickness.push_back(Thickness);
+ m_AC_Material.push_back(Material);
+ m_AC_Angle.push_back(Angle);
+}
+
+void ZDD::Add_Entry_Exit(G4double Z, double Thickness, string Material) {
+ m_Entry_Exit_Z.push_back(Z);
+ m_Entry_Exit_Thickness.push_back(Thickness);
+ m_Entry_Exit_Material.push_back(Material);
+}
+
+void ZDD::Add_Plastic(string Material, G4double Width, double Length,
+ double Thickness, G4ThreeVector Pos) {
+ m_Plastic_Material.push_back(Material);
+ m_Plastic_Width.push_back(Width);
+ m_Plastic_Length.push_back(Length);
+ m_Plastic_Thickness.push_back(Thickness);
+ m_Plastic_Position.push_back(Pos);
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+/*G4LogicalVolume* ZDD::BuildSquareDetector(){
+ if(!m_SquareDetector){
+ G4Box* box = new G4Box("ZDD_Box",ZDD_NS::Width*0.5,
+ ZDD_NS::Width*0.5,ZDD_NS::Thickness*0.5);
+
+ G4Material* DetectorMaterial = MaterialManager::getInstance()->GetMaterialFromLibrary(ZDD_NS::Material);
+ m_SquareDetector = new G4LogicalVolume(box,DetectorMaterial,"logic_ZDD_Box",0,0,0);
+ m_SquareDetector->SetVisAttributes(m_VisSquare);
+ m_SquareDetector->SetSensitiveDetector(m_ZDDScorer);
+ }
+ return m_SquareDetector;
+}*/
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+/*G4LogicalVolume* ZDD::BuildCylindricalDetector(){
+ if(!m_CylindricalDetector){
+ G4Tubs* tub = new G4Tubs("ZDD_Cyl",0,ZDD_NS::Radius,ZDD_NS::Thickness*0.5,0,360*deg);
+
+ G4Material* DetectorMaterial = MaterialManager::getInstance()->GetMaterialFromLibrary(ZDD_NS::Material);
+ m_CylindricalDetector = new G4LogicalVolume(tub,DetectorMaterial,"logic_ZDD_tub",0,0,0);
+ m_CylindricalDetector->SetVisAttributes(m_VisSquare);
+ m_CylindricalDetector->SetSensitiveDetector(m_ZDDScorer);
+
+ }
+ return m_CylindricalDetector;
+}*/
+
+G4LogicalVolume* ZDD::Build_Drift_Chamber_1() {
+ if (!m_Drift_Chamber_1) {
+ G4Box* box = new G4Box("ZDD_Drift_Chamber_1", Drift_Chamber_Width * 0.5, Drift_Chamber_Length * 0.5,
+ m_Drift_Chamber_Thickness[0] * 0.5);
+
+ G4Material* DetectorMaterial
+ = MaterialManager::getInstance()->GetGasFromLibrary(
+ m_Drift_Chamber_Gas[0], m_Drift_Chamber_Pressure[0], m_Drift_Chamber_Temperature[0]);
+
+ m_Drift_Chamber_1 = new G4LogicalVolume(box, DetectorMaterial, "logic_ZDD_Drift_Chamber_1", 0, 0,
+ 0);
+ m_Drift_Chamber_1->SetVisAttributes(m_Vis_Drift_Chamber_Gas);
+ m_Drift_Chamber_1->SetSensitiveDetector(m_Drift_Chamber_Scorer_1);
+ }
+ return m_Drift_Chamber_1;
+
+}
+
+G4LogicalVolume* ZDD::Build_Drift_Chamber_2() {
+ if (!m_Drift_Chamber_2) {
+ G4Box* box = new G4Box("ZDD_Drift_Chamber_2", Drift_Chamber_Width * 0.5, Drift_Chamber_Length * 0.5,
+ m_Drift_Chamber_Thickness[1] * 0.5);
+
+ G4Material* DetectorMaterial
+ = MaterialManager::getInstance()->GetGasFromLibrary(
+ m_Drift_Chamber_Gas[1], m_Drift_Chamber_Pressure[1], m_Drift_Chamber_Temperature[1]);
+
+ m_Drift_Chamber_2 = new G4LogicalVolume(box, DetectorMaterial, "logic_ZDD_Drift_Chamber_2", 0, 0,
+ 0);
+ m_Drift_Chamber_2->SetVisAttributes(m_Vis_Drift_Chamber_Gas);
+ m_Drift_Chamber_2->SetSensitiveDetector(m_Drift_Chamber_Scorer_2);
+ }
+ return m_Drift_Chamber_2;
+}
+
+void ZDD::ClearGeometry(){
+
+ m_Drift_Chamber_Z.clear();
+ m_Drift_Chamber_Gas.clear();
+ m_Drift_Chamber_Pressure.clear();
+ m_Drift_Chamber_Temperature.clear();
+ m_Drift_Chamber_Thickness.clear();
+
+ m_Ionisation_Chamber_Z.clear();
+ m_Ionisation_Chamber_Gas.clear();
+ m_Ionisation_Chamber_Pressure.clear();
+ m_Ionisation_Chamber_Temperature.clear();
+ m_Ionisation_Chamber_Thickness.clear();
+
+ m_Gas_Gap_Z.clear();
+ m_Gas_Gap_Gas.clear();
+ m_Gas_Gap_Pressure.clear();
+ m_Gas_Gap_Temperature.clear();
+ m_Gas_Gap_Thickness.clear();
+
+ m_AC_Material.clear();
+ m_AC_Rotation.clear();
+ m_AC_Thickness.clear();
+ m_AC_Z.clear();
+ m_AC_Angle.clear();
+
+ m_Entry_Exit_Material.clear();
+ m_Entry_Exit_Z.clear();
+ m_Entry_Exit_Thickness.clear();
+
+ m_Plastic_Length.clear();
+ m_Plastic_Width.clear();
+ m_Plastic_Thickness.clear();
+ m_Plastic_Position.clear();
+ m_Plastic_Material.clear();
+}
+//....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 ZDD::ReadConfiguration(NPL::InputParser parser){
+ //ClearGeometry();
+
+ vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("ZDD");
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << "//// " << blocks.size() << " detectors found " << endl;
+
+ vector<string> TokenZDD = {"R", "Theta"};
+ vector<string> TokenDC = {"Z","Thickness", "Gas", "Pressure", "Temperature"};
+ vector<string> TokenIC = {"Z", "Thickness", "Gas", "Pressure", "Temperature"};
+ vector<string> TokenAC = {"Z", "Thickness", "Material"};
+ vector<string> TokenEntryExit = {"Z", "Thickness", "Material"};
+ vector<string> TokenPlastic = {"Material", "Width", "Length", "Thickness", "Pos"};
+
+ for(unsigned int i = 0 ; i < blocks.size() ; i++){
+ /*if(blocks[i]->HasTokenList(cart)){
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// ZDD " << i+1 << endl;
+
+ G4ThreeVector Pos = NPS::ConvertVector(blocks[i]->GetTVector3("POS","mm"));
+ string Shape = blocks[i]->GetString("Shape");
+ AddDetector(Pos,Shape);
+ }
+ else if(blocks[i]->HasTokenList(sphe)){
+ if(NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// ZDD " << i+1 << endl;
+ double R = blocks[i]->GetDouble("R","mm");
+ double Theta = blocks[i]->GetDouble("Theta","deg");
+ double Phi = blocks[i]->GetDouble("Phi","deg");
+ string Shape = blocks[i]->GetString("Shape");
+ AddDetector(R,Theta,Phi,Shape);
+ }*/
+ if (blocks[i]->HasTokenList(TokenZDD)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// ZDD " << i + 1 << endl;
+ G4double R = blocks[i]->GetDouble("R", "mm");
+ G4double Theta = blocks[i]->GetDouble("Theta", "deg");
+ Add_ZDD(R, Theta);
+ }
+ else if (blocks[i]->GetMainValue() == "DC"
+ && blocks[i]->HasTokenList(TokenDC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// DC" << i + 1 << endl;
+ G4double Z = blocks[i]->GetDouble("Z", "mm");
+ G4double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ string Gas = blocks[i]->GetString("Gas");
+ G4double Pressure = blocks[i]->GetDouble("Pressure", "bar");
+ G4double Temperature = blocks[i]->GetDouble("Temperature", "kelvin");
+ Add_Drift_Chamber(Z, Thickness, Gas, Pressure, Temperature);
+ }
+ else if (blocks[i]->GetMainValue() == "IC"
+ && blocks[i]->HasTokenList(TokenIC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// IC" << ICcounter+1 << endl;
+ G4double Z = blocks[i]->GetDouble("Z", "mm");
+ G4double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ string Gas = blocks[i]->GetString("Gas");
+ G4double Pressure = blocks[i]->GetDouble("Pressure", "bar");
+ G4double Temperature = blocks[i]->GetDouble("Temperature", "kelvin");
+ Add_Ionisation_Chamber(Z, Thickness, Gas, Pressure, Temperature);
+ ICcounter++;
+ }
+ else if (blocks[i]->GetMainValue() == "GasGap"
+ && blocks[i]->HasTokenList(TokenIC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// GasGap" << GGcounter+1 << endl;
+ G4double Z = blocks[i]->GetDouble("Z", "mm");
+ G4double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ string Gas = blocks[i]->GetString("Gas");
+ G4double Pressure = blocks[i]->GetDouble("Pressure", "bar");
+ G4double Temperature = blocks[i]->GetDouble("Temperature", "kelvin");
+ Add_Gas_Gap(Z, Thickness, Gas, Pressure, Temperature);
+ GGcounter++;
+ }
+ else if (blocks[i]->GetMainValue() == "AC"
+ && blocks[i]->HasTokenList(TokenAC)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// AC" << ACcounter+1 << endl;
+ G4double Z = blocks[i]->GetDouble("Z", "mm");
+ G4double Thickness = blocks[i]->GetDouble("Thickness", "um");
+ string Material = blocks[i]->GetString("Material");
+ G4double Angle = blocks[i]->GetDouble("Theta","deg");
+ Add_AC(Z, Thickness, Material, Angle);
+ ACcounter++;
+ }
+ else if (blocks[i]->GetMainValue() == "EntryExit"
+ && blocks[i]->HasTokenList(TokenEntryExit)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// AC" << Entry_Exit_counter+1 << endl;
+ G4double Z = blocks[i]->GetDouble("Z", "mm");
+ G4double Thickness = blocks[i]->GetDouble("Thickness", "um");
+ string Material = blocks[i]->GetString("Material");
+ Add_Entry_Exit(Z, Thickness, Material);
+ Entry_Exit_counter++;
+ }
+ else if (blocks[i]->GetMainValue() == "Plastic"
+ && blocks[i]->HasTokenList(TokenPlastic)) {
+ if (NPOptionManager::getInstance()->GetVerboseLevel())
+ cout << endl << "//// Plastic" << Plasticcounter + 1 << endl;
+ string Material = blocks[i]->GetString("Material");
+ G4double Width = blocks[i]->GetDouble("Width", "mm");
+ G4double Length = blocks[i]->GetDouble("Length", "mm");
+ G4double Thickness = blocks[i]->GetDouble("Thickness", "mm");
+ G4ThreeVector Pos = NPS::ConvertVector(blocks[i]->GetTVector3("Pos", "mm"));
+ Add_Plastic(Material, Width, Length, Thickness, Pos);
+ Plasticcounter++;
+ }
+ 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 ZDD::ConstructDetector(G4LogicalVolume* world){
+// for (unsigned short i = 0 ; i < m_R.size() ; i++) {
+
+ // Mother Volume of ZDD
+ G4double R = m_R + ZDD_Volume_Thickness* 0.5;
+
+ G4double X = R * sin(m_Theta) * cos(Phi);
+ G4double Y = R * sin(m_Theta) * sin(Phi);
+ G4double Z = R * cos(m_Theta);
+ G4ThreeVector Det_pos = G4ThreeVector(X, Y, Z);
+
+ // Motehr Volume de la ZDD
+ G4RotationMatrix* Rot1 = new G4RotationMatrix();
+ Rot1->rotateY(m_Theta);
+
+ G4Box* MotherSolid
+ = new G4Box("MotherVolume", ZDD_Volume_Width * 0.5,
+ ZDD_Volume_Length * 0.5, ZDD_Volume_Thickness* 0.5);
+
+ G4Material* VolumeMaterial
+ = MaterialManager::getInstance()->GetMaterialFromLibrary("Vacuum");
+ G4LogicalVolume* MotherVolume = new G4LogicalVolume(
+ MotherSolid, VolumeMaterial, "MotherVolume", 0, 0, 0);
+
+ new G4PVPlacement(G4Transform3D(*Rot1, Det_pos), MotherVolume, "MotherVolume",
+ world, false, 0);
+ MotherVolume->SetVisAttributes(m_Vis_ZDD);
+
+
+ // Drift Chamber 1
+ if(m_Drift_Chamber_Z.size() > 0){
+ new G4PVPlacement(0, G4ThreeVector(0, 0, (-ZDD_Volume_Thickness + m_Drift_Chamber_Thickness[0]) * 0.5 + m_Drift_Chamber_Z[0]),
+ Build_Drift_Chamber_1(), "Drift_Chamber_1", MotherVolume, false, 1);
+
+ // Drift Chamber 2
+ new G4PVPlacement(0, G4ThreeVector(0, 0, (-ZDD_Volume_Thickness + m_Drift_Chamber_Thickness[1])* 0.5 + m_Drift_Chamber_Z[1]),
+ Build_Drift_Chamber_2(), "Drift_Chamber_2", MotherVolume, false, 2);
+ }
+
+ // Ionisation Chambers (gas)
+ G4RotationMatrix* Rot_Ionisation_Chamber = new G4RotationMatrix();
+ Rot_Ionisation_Chamber->rotateX(Ionisation_Chamber_Angle);
+
+ for (int i = 0; i < ICcounter; i++) {
+ G4Box* box = new G4Box("ZDD_IC", Ionisation_Chamber_Width * 0.5, Ionisation_Chamber_Length * 0.5,
+ m_Ionisation_Chamber_Thickness[i] * 0.5);
+ G4Material* GasIC = MaterialManager::getInstance()->GetGasFromLibrary(
+ m_Ionisation_Chamber_Gas[i], m_Ionisation_Chamber_Pressure[i], m_Ionisation_Chamber_Temperature[i]);
+ G4LogicalVolume* IC
+ = new G4LogicalVolume(box, GasIC, "logic_ZDD_IC", 0, 0, 0);
+ IC->SetVisAttributes(m_Vis_Ionisation_Chamber_Gas);
+ IC->SetSensitiveDetector(m_IC_Scorer);
+ new G4PVPlacement(G4Transform3D(*Rot_Ionisation_Chamber,
+ G4ThreeVector(0, 0, (-ZDD_Volume_Thickness + m_Ionisation_Chamber_Thickness[i])*0.5 + m_Ionisation_Chamber_Z[i])),
+ IC, "IC", MotherVolume, false, i);
+ }
+ // gas Gap
+ for (int i = 0; i < GGcounter; i++) {
+ G4Box* box = new G4Box("ZDD_GG", Gas_Gap_Width * 0.5, Gas_Gap_Length * 0.5,
+ m_Gas_Gap_Thickness[i] * 0.5);
+ G4Material* GasGG = MaterialManager::getInstance()->GetGasFromLibrary(
+ m_Gas_Gap_Gas[i], m_Gas_Gap_Pressure[i], m_Gas_Gap_Temperature[i]);
+ G4LogicalVolume* GG
+ = new G4LogicalVolume(box, GasGG, "logic_ZDD_GG", 0, 0, 0);
+ GG->SetVisAttributes(m_Vis_Ionisation_Chamber_Gas);
+ new G4PVPlacement(G4Transform3D(*Rot_Ionisation_Chamber,
+ G4ThreeVector(0, 0, (-ZDD_Volume_Thickness + m_Gas_Gap_Thickness[i])*0.5 + m_Gas_Gap_Z[i])),
+ GG, "GG", MotherVolume, false, i);
+ }
+
+ // AC
+
+ for (int i = 0; i < ACcounter; i++) {
+ G4RotationMatrix* Rot_AC = new G4RotationMatrix();
+ Rot_AC->rotateX(m_AC_Angle[i]);
+ G4Box* box = new G4Box("ZDD_AC", AC_Width * 0.5, AC_Length * 0.5,
+ m_AC_Thickness[i] * 0.5);
+ G4Material* MaterialAC = MaterialManager::getInstance()->GetMaterialFromLibrary(m_AC_Material[i]);
+ G4LogicalVolume* AC
+ = new G4LogicalVolume(box, MaterialAC, "logic_ZDD_AC", 0, 0, 0);
+ AC->SetVisAttributes(m_Vis_AC);
+ new G4PVPlacement(G4Transform3D(*Rot_AC,
+ G4ThreeVector(0, 0, (-ZDD_Volume_Thickness + m_AC_Thickness[i])*0.5 + m_AC_Z[i])),
+ AC, "AC", MotherVolume, false, i);
+ }
+
+ // Entry_Exit
+ for (int i = 0; i < Entry_Exit_counter; i++) {
+ G4Box* box = new G4Box("ZDD_Entry_Exit", Kapton_Foil_Width * 0.5, Kapton_Foil_Length * 0.5,
+ m_Entry_Exit_Thickness[i] * 0.5);
+ G4Material* Material_Entry_Exit = MaterialManager::getInstance()->GetMaterialFromLibrary(m_Entry_Exit_Material[i]);
+ G4LogicalVolume* Entry_Exit
+ = new G4LogicalVolume(box, Material_Entry_Exit, "logic_ZDD_Entry_Exit", 0, 0, 0);
+ Entry_Exit->SetVisAttributes(m_Vis_AC);
+ new G4PVPlacement(G4Transform3D(*Rot_Ionisation_Chamber,
+ G4ThreeVector(0, 0, (-ZDD_Volume_Thickness + m_Entry_Exit_Thickness[i])*0.5 + m_Entry_Exit_Z[i])),
+ Entry_Exit, "Entry_Exit", MotherVolume, false, i);
+ }
+
+ // Plastic
+ for (int i = 0; i < Plasticcounter; i++) {
+ G4Box* box = new G4Box("ZDD_Plastic", m_Plastic_Width[i] * 0.5, m_Plastic_Length[i] * 0.5,
+ m_Plastic_Thickness[i] * 0.5);
+ G4Material* Material = MaterialManager::getInstance()->GetMaterialFromLibrary(m_Plastic_Material[i]);
+ G4LogicalVolume* Plastic
+ = new G4LogicalVolume(box, Material, "logic_ZDD_Plastic", 0, 0, 0);
+ Plastic->SetVisAttributes(m_Vis_Plastic);
+ Plastic->SetSensitiveDetector(m_Plastic_Scorer);
+ new G4PVPlacement(
+ 0,
+ G4ThreeVector(m_Plastic_Position[i][0], m_Plastic_Position[i][1],
+ (-ZDD_Volume_Thickness + m_Plastic_Thickness[i])*0.5 + m_Plastic_Position[i][2]),
+ Plastic, "Plastic", MotherVolume, false, 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()*ZDD_NS::Thickness*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);
+
+ if(m_Shape[i] == "Cylindrical"){
+ new G4PVPlacement(G4Transform3D(*Rot,Det_pos),
+ BuildCylindricalDetector(),
+ "ZDD",world,false,i+1);
+ }
+
+ else if(m_Shape[i] == "Square"){
+ new G4PVPlacement(G4Transform3D(*Rot,Det_pos),
+ BuildSquareDetector(),
+ "ZDD",world,false,i+1);
+ }*/
+ // }
+}
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+// Add Detector branch to the EventTree.
+// Called After DetecorConstruction::AddDetector Method
+void ZDD::InitializeRootOutput(){
+ RootOutput *pAnalysis = RootOutput::getInstance();
+ TTree *pTree = pAnalysis->GetTree();
+ if(!pTree->FindBranch("ZDD")){
+ pTree->Branch("ZDD", "TZDDData", &m_Event) ;
+ }
+ pTree->SetBranchAddress("ZDD", &m_Event) ;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+// Read sensitive part and fill the Root tree.
+// Called at in the EventAction::EndOfEventAvtion
+void ZDD::ReadSensitive(const G4Event* ){
+ m_Event->Clear();
+ ///////////
+ // Calorimeter scorer
+ CalorimeterScorers::PS_Calorimeter* Scorer_IC
+ = (CalorimeterScorers::PS_Calorimeter*)m_IC_Scorer->GetPrimitive(0);
+ unsigned int size_IC = Scorer_IC->GetMult();
+ for (unsigned int i = 0; i < size_IC; i++) {
+ vector<unsigned int> level = Scorer_IC->GetLevel(i);
+ G4double Energy
+ = RandGauss::shoot(Scorer_IC->GetEnergy(i), Scorer_IC->GetEnergy(i)*ZDD_NS::ResoEnergy_IC/(100*2.355));
+ if (Energy > ZDD_NS::EnergyThreshold_IC) {
+ G4double Time = RandGauss::shoot(Scorer_IC->GetTime(i), ZDD_NS::ResoTime_IC);
+ int DetectorNbr = level[0];
+ m_Event->Set_IC_Energy(DetectorNbr, Energy);
+ m_Event->Set_IC_Time(DetectorNbr, Time);
+ }
+ }
+
+ CalorimeterScorers::PS_Calorimeter* Scorer_Plastic
+ = (CalorimeterScorers::PS_Calorimeter*)m_Plastic_Scorer->GetPrimitive(
+ 0);
+ unsigned int size_Plastic = Scorer_Plastic->GetMult();
+ for (unsigned int i = 0; i < size_Plastic; i++) {
+ vector<unsigned int> level = Scorer_Plastic->GetLevel(i);
+ G4double Energy
+ = RandGauss::shoot(Scorer_Plastic->GetEnergy(i), Scorer_Plastic->GetEnergy(i)*ZDD_NS::ResoEnergy_Plastic/(100*2.355));
+ if (Energy > ZDD_NS::EnergyThreshold_Plastic) {
+ G4double Time = RandGauss::shoot(Scorer_Plastic->GetTime(i), ZDD_NS::ResoTime_Plastic);
+ int DetectorNbr = level[0];
+ m_Event->Set_Plastic_Energy(DetectorNbr, Energy);
+ m_Event->Set_Plastic_Time(DetectorNbr, Time);
+ }
+ }
+
+ ///////////
+ // Calorimeter scorer
+ /*CalorimeterScorers::PS_Calorimeter* Scorer= (CalorimeterScorers::PS_Calorimeter*) m_ZDDScorer->GetPrimitive(0);
+
+ unsigned int size = Scorer->GetMult();
+ for(unsigned int i = 0 ; i < size ; i++){
+ vector<unsigned int> level = Scorer->GetLevel(i);
+ double Energy = RandGauss::shoot(Scorer->GetEnergy(i),ZDD_NS::ResoEnergy);
+ if(Energy>ZDD_NS::EnergyThreshold){
+ double Time = RandGauss::shoot(Scorer->GetTime(i),ZDD_NS::ResoTime);
+ int DetectorNbr = level[0];
+ m_Event->SetEnergy(DetectorNbr,Energy);
+ m_Event->SetTime(DetectorNbr,Time);
+ }
+ }*/
+ DriftChamberScorers::PS_DriftChamber* Scorer_DC_1
+ = (DriftChamberScorers::PS_DriftChamber*)m_Drift_Chamber_Scorer_1->GetPrimitive(0);
+ unsigned int size1 = Scorer_DC_1->GetMult();
+ for (unsigned int i = 0; i < size1; i++) {
+ vector<unsigned int> level = Scorer_DC_1->GetLevel(i);
+ G4double DriftTime
+ = RandGauss::shoot(Scorer_DC_1->GetDriftTime(i)/Scorer_DC_1->GetCounter(i), ZDD_NS::ResoDriftTime);
+ int DetectorNbr = level[0];
+ m_Event->Set_DC_Time(DetectorNbr, DriftTime);
+ }
+
+ DriftChamberScorers::PS_DriftChamber* Scorer_DC_2
+ = (DriftChamberScorers::PS_DriftChamber*)m_Drift_Chamber_Scorer_2->GetPrimitive(0);
+ unsigned int size2 = Scorer_DC_2->GetMult();
+ for (unsigned int i = 0; i < size2; i++) {
+ vector<unsigned int> level = Scorer_DC_2->GetLevel(i);
+ G4double DriftTime
+ = RandGauss::shoot(Scorer_DC_2->GetDriftTime(i)/Scorer_DC_2->GetCounter(i), ZDD_NS::ResoDriftTime);
+ int DetectorNbr = level[0];
+ m_Event->Set_DC_Time(DetectorNbr, DriftTime);
+ }
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+////////////////////////////////////////////////////////////////
+void ZDD::InitializeScorers() {
+ // This check is necessary in case the geometry is reloaded
+ bool already_exist = false;
+ m_Drift_Chamber_Scorer_1 = CheckScorer("Drift_Chamber_Scorer_1", already_exist);
+ m_Drift_Chamber_Scorer_2 = CheckScorer("Drift_Chamber_Scorer_2", already_exist);
+ m_IC_Scorer = CheckScorer("IC_Scorer",already_exist) ;
+ m_Plastic_Scorer = CheckScorer("Plastic_Scorer",already_exist) ;
+
+ if(already_exist)
+ return ;
+
+ // Otherwise the scorer is initialised
+ vector<int> level;
+ level.push_back(0);
+
+ //G4VPrimitiveScorer* Calorimeter= new CalorimeterScorers::PS_Calorimeter("Calorimeter",level, 0) ;
+ //G4VPrimitiveScorer* Interaction= new InteractionScorers::PS_Interactions("Interaction",ms_InterCoord, 0) ;
+ G4VPrimitiveScorer* Drift1 = new DriftChamberScorers::PS_DriftChamber("Drift", level, DriftDir1, DriftSpeed, 0);
+ G4VPrimitiveScorer* Drift2 = new DriftChamberScorers::PS_DriftChamber("Drift", level, DriftDir2, DriftSpeed, 0);
+ G4VPrimitiveScorer* IC = new CalorimeterScorers::PS_Calorimeter("IC", level, 0);
+ G4VPrimitiveScorer* Plastic = new CalorimeterScorers::PS_Calorimeter("Plastic", level, 0);
+
+
+ //and register it to the multifunctionnal detector
+
+ //m_ZDDScorer->RegisterPrimitive(Calorimeter);
+ //m_ZDDScorer->RegisterPrimitive(Interaction);
+ m_Drift_Chamber_Scorer_1->RegisterPrimitive(Drift1);
+ m_Drift_Chamber_Scorer_2->RegisterPrimitive(Drift2);
+ m_IC_Scorer->RegisterPrimitive(IC);
+ m_Plastic_Scorer->RegisterPrimitive(Plastic);
+
+ G4SDManager::GetSDMpointer()->AddNewDetector(m_Drift_Chamber_Scorer_1) ;
+ G4SDManager::GetSDMpointer()->AddNewDetector(m_Drift_Chamber_Scorer_2) ;
+ G4SDManager::GetSDMpointer()->AddNewDetector(m_IC_Scorer);
+ G4SDManager::GetSDMpointer()->AddNewDetector(m_Plastic_Scorer);
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+////////////////////////////////////////////////////////////////////////////////
+// Construct Method to be pass to the DetectorFactory //
+////////////////////////////////////////////////////////////////////////////////
+NPS::VDetector* ZDD::Construct(){
+ return (NPS::VDetector*) new ZDD();
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+////////////////////////////////////////////////////////////////////////////////
+// Registering the construct method to the factory //
+////////////////////////////////////////////////////////////////////////////////
+extern"C" {
+ class proxy_nps_ZDD{
+ public:
+ proxy_nps_ZDD(){
+ NPS::DetectorFactory::getInstance()->AddToken("ZDD","ZDD");
+ NPS::DetectorFactory::getInstance()->AddDetector("ZDD",ZDD::Construct);
+ }
+ };
+
+ proxy_nps_ZDD p_nps_ZDD;
+}
diff --git a/NPSimulation/Detectors/ZDD/ZDD.hh b/NPSimulation/Detectors/ZDD/ZDD.hh
new file mode 100644
index 0000000000000000000000000000000000000000..7389a15ed19eb29ba3b04ed1e7ba1c439503fac9
--- /dev/null
+++ b/NPSimulation/Detectors/ZDD/ZDD.hh
@@ -0,0 +1,196 @@
+#ifndef ZDD_h
+#define ZDD_h 1
+/*****************************************************************************
+ * Copyright (C) 2009-2022 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: Hugo Jacob contact address: hjacob@ijclab.in2p3.fr *
+ * *
+ * Creation Date : October 2022 *
+ * Last update : *
+ *---------------------------------------------------------------------------*
+ * Decription: *
+ * This class describe ZDD 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 "TZDDData.h"
+#include "NPInputParser.h"
+
+class ZDD : public NPS::VDetector{
+ ////////////////////////////////////////////////////
+ /////// Default Constructor and Destructor /////////
+ ////////////////////////////////////////////////////
+ public:
+ ZDD() ;
+ virtual ~ZDD() ;
+
+ ////////////////////////////////////////////////////
+ /////// Specific Function of this Class ///////////
+ ////////////////////////////////////////////////////
+ public:
+ // Cartesian
+ void AddDetector(G4ThreeVector POS, string Shape);
+ // Spherical
+ void AddDetector(double R,double Theta,double Phi,string Shape);
+
+ void Add_Drift_Chamber(G4double Z,double Thickness, string Gas, double Pressure, double Temperature);
+
+ void Add_Ionisation_Chamber(G4double Z, double Thickness, string Gas, double Pressure, double Temperature);
+
+ void Add_Gas_Gap(G4double Z, double Thickness, string Gas, double Pressure, double Temperature);
+
+ void Add_AC(G4double Z, double Thickness, string Material, G4double Angle);
+
+ void Add_Entry_Exit(G4double Z, double Thickness, string Material);
+
+ void Add_Plastic(string Material, G4double Width, double Length, double Thickness, G4ThreeVector Pos);
+
+ void Add_ZDD(G4double R, double theta);
+
+
+ G4LogicalVolume* BuildSquareDetector();
+ G4LogicalVolume* BuildCylindricalDetector();
+ G4LogicalVolume* Build_Drift_Chamber_1();
+ G4LogicalVolume* Build_Drift_Chamber_2();
+
+ private:
+ G4LogicalVolume* m_Drift_Chamber_1;
+ G4LogicalVolume* m_Drift_Chamber_2;
+
+ G4double ICcounter;
+ G4double ACcounter;
+ G4double GGcounter;
+ G4double Entry_Exit_counter;
+ G4double Plasticcounter;
+ //G4LogicalVolume* m_SquareDetector;
+ //G4LogicalVolume* m_CylindricalDetector;
+
+ ////////////////////////////////////////////////////
+ ////// Inherite from NPS::VDetector class /////////
+ ////////////////////////////////////////////////////
+ public:
+
+ void ClearGeometry();
+
+ // 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_ZDDScorer ;
+ G4MultiFunctionalDetector* m_Drift_Chamber_Scorer_2 ;
+ G4MultiFunctionalDetector* m_Drift_Chamber_Scorer_1 ;
+ G4MultiFunctionalDetector* m_IC_Scorer ;
+ G4MultiFunctionalDetector* m_Plastic_Scorer ;
+
+ ////////////////////////////////////////////////////
+ ///////////Event class to store Data////////////////
+ ////////////////////////////////////////////////////
+ private:
+ TZDDData* m_Event ;
+
+ ////////////////////////////////////////////////////
+ ///////////////Private intern Data//////////////////
+ ////////////////////////////////////////////////////
+ private: // Geometry
+ // Detector Coordinate
+ G4double m_R;
+ G4double m_Theta;
+ //vector<double> m_Phi;
+
+ // Shape type
+ vector<string> m_Shape ;
+
+ // Visualisation Attribute
+ G4VisAttributes* m_VisSquare;
+ G4VisAttributes* m_VisCylinder;
+
+ // Drift Chambers Gas, pressure etc.
+ vector<G4double> m_Drift_Chamber_Z;
+ vector<string> m_Drift_Chamber_Gas;
+ vector<G4double> m_Drift_Chamber_Pressure;
+ vector<G4double> m_Drift_Chamber_Temperature;
+ vector<G4double> m_Drift_Chamber_Thickness;
+
+ // Ionisation Chambers Gas, pressure etc.
+ vector<G4double> m_Ionisation_Chamber_Z;
+ vector<G4double> m_Ionisation_Chamber_Thickness;
+ vector<string> m_Ionisation_Chamber_Gas;
+ vector<G4double> m_Ionisation_Chamber_Pressure;
+ vector<G4double> m_Ionisation_Chamber_Temperature;
+
+ // Gas Gap
+ vector<G4double> m_Gas_Gap_Z;
+ vector<G4double> m_Gas_Gap_Thickness;
+ vector<string> m_Gas_Gap_Gas;
+ vector<G4double> m_Gas_Gap_Pressure;
+ vector<G4double> m_Gas_Gap_Temperature;
+
+ // Anode Cathodes
+ vector<G4double> m_AC_Z;
+ vector<G4bool> m_AC_Rotation;
+ vector<G4double> m_AC_Thickness;
+ vector<string> m_AC_Material;
+ vector<G4double> m_AC_Angle;
+
+ // Entry and Exit (Kapton)
+ vector<G4double> m_Entry_Exit_Z;
+ vector<G4double> m_Entry_Exit_Thickness;
+ vector<string> m_Entry_Exit_Material;
+
+ // Plastic Pos, size etc.
+ vector<G4double> m_Plastic_Width;
+ vector<G4double> m_Plastic_Length;
+ vector<G4double> m_Plastic_Thickness;
+ vector<string> m_Plastic_Material;
+ vector<G4ThreeVector> m_Plastic_Position;
+
+
+ // Visualisation Attributes;
+ G4VisAttributes* m_Vis_Drift_Chamber_Gas;
+ G4VisAttributes* m_Vis_Ionisation_Chamber_Gas;
+ G4VisAttributes* m_Vis_AC;
+ G4VisAttributes* m_Vis_Plastic;
+ G4VisAttributes* m_Vis_ZDD;
+ // Needed for dynamic loading of the library
+ public:
+ static NPS::VDetector* Construct();
+};
+#endif
diff --git a/Projects/MUST_AND_ZDD/Analysis.cxx b/Projects/MUST_AND_ZDD/Analysis.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..92feddcf92203f93f181f24886655dd1af5a2dcc
--- /dev/null
+++ b/Projects/MUST_AND_ZDD/Analysis.cxx
@@ -0,0 +1,433 @@
+/*****************************************************************************
+ * 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 MUST_AND_ZDD analysis project *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ *****************************************************************************/
+
+#include<iostream>
+using namespace std;
+#include"Analysis.h"
+#include"NPAnalysisFactory.h"
+#include"NPDetectorManager.h"
+////////////////////////////////////////////////////////////////////////////////
+Analysis::Analysis(){
+}
+////////////////////////////////////////////////////////////////////////////////
+Analysis::~Analysis(){
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::Init(){
+ InitialConditions = new TInitialConditions;
+ ReactionConditions = new TReactionConditions;
+ InitOutputBranch();
+ InitInputBranch();
+ ZDD= (TZDDPhysics*) m_DetectorManager->GetDetector("ZDD");
+ M2 = (TMust2Physics*) m_DetectorManager -> GetDetector("M2Telescope");
+ reaction->ReadConfigurationFile(NPOptionManager::getInstance()->GetReactionFile());
+ OriginalBeamEnergy = reaction->GetBeamEnergy();
+
+
+ string Path = "../../Inputs/EnergyLoss/";
+ string TargetMaterial = m_DetectorManager->GetTargetMaterial();
+ TargetThickness = m_DetectorManager->GetTargetThickness();
+ string beam= NPL::ChangeNameToG4Standard(reaction->GetNucleus1()->GetName());
+ string heavy_ejectile= NPL::ChangeNameToG4Standard(reaction->GetNucleus4()->GetName());
+ string light=NPL::ChangeNameToG4Standard(reaction->GetNucleus3()->GetName());
+
+
+ Beam_Target = NPL::EnergyLoss(beam+"_"+TargetMaterial+".G4table","G4Table",100);
+ Heavy_Target = NPL::EnergyLoss(heavy_ejectile+"_"+TargetMaterial+".G4table","G4Table",100);
+ LightTarget = NPL::EnergyLoss(light+"_"+TargetMaterial+".G4table","G4Table",100 );
+ LightAl = NPL::EnergyLoss(light+"_Al.G4table","G4Table",100);
+ LightSi = NPL::EnergyLoss(light+"_Si.G4table","G4Table",100);
+
+ ///////////////////////// Creating EnerlyLoss Objects for ZDD ////////////////////////////////////
+ for(int i = 0; i < ZDD->Get_AC_Material().size(); i++){
+ Heavy_IC_Mylar.push_back(NPL::EnergyLoss(heavy_ejectile+"_"+ZDD->Get_AC_Material()[i]+".G4table","G4Table",100));
+ }
+ for(int i = 0; i < ZDD->Get_Entry_Exit_Material().size(); i++){
+ Heavy_IC_Windows.push_back(NPL::EnergyLoss(heavy_ejectile+"_"+ZDD->Get_Entry_Exit_Material()[i]+".G4table","G4Table",100));
+ }
+ for(int i = 0; i < ZDD->Get_m_Gas_Gap_Gas().size(); i++){
+ double pressure = ZDD->Get_m_Gas_Gap_Pressure()[i]/bar;
+ string pres = to_string(pressure);
+ pres.erase(pres.find_last_not_of('0')+1,pres.size());
+ if(pressure >= 1){
+ pres.erase(pres.find_last_not_of('.')+1,pres.size());
+ }
+ int Temp = ZDD->Get_m_Gas_Gap_Temperature()[i];
+ Heavy_IC_Gas.push_back(NPL::EnergyLoss(heavy_ejectile+"_"+ZDD->Get_m_Gas_Gap_Gas()[i]+"_"+pres+"bar_"+to_string(Temp)+"K.G4table","G4Table",100));
+
+ }
+ for(int i = 0; i < ZDD->Get_m_Drift_Chamber_Gas().size(); i++){
+ double pressure = ZDD->Get_m_Drift_Chamber_Pressure()[i]/bar;
+ string pres = to_string(pressure);
+ pres.erase(pres.find_last_not_of('0')+1,pres.size());
+ if(pressure >= 1){
+ pres.erase(pres.find_last_not_of('.')+1,pres.size());
+ }
+ int Temp = ZDD->Get_m_Drift_Chamber_Temperature()[i];
+ Heavy_DC_Gas.push_back(NPL::EnergyLoss(heavy_ejectile+"_"+ZDD->Get_m_Drift_Chamber_Gas()[i]+"_"+pres+"bar_"+to_string(Temp)+"K.G4table","G4Table",100));
+
+ }
+ ///////////////////////// Creating EnergyLoss Objects for MUST2 ////////////////////////////////////
+
+ //////////////////////////// Treating Drift Chamber (for this part, we assume the drift speed, in further analysis we shall use CATS)
+ //////////////////////////// (information to retrieve drift coefficients)
+
+
+ ///////////////////////////// Initialize some important parameters //////////////////////////////////
+
+ Drift_Speed = 1 * cm / microsecond;
+ ZDir.SetXYZ(0.,0.,1.);
+
+
+ DetectorNumber = 0;
+ ThetaNormalTarget = 0;
+ ThetaM2Surface = 0;
+ ThetaMGSurface = 0;
+ Si_E_M2 = 0;
+ CsI_E_M2 = 0;
+ Energy = 0;
+ ThetaGDSurface = 0;
+ M2_X = 0;
+ M2_Y = 0;
+ M2_Z = 0;
+ M2_dE = 0;
+ BeamDirection = TVector3(0,0,1);
+}
+
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::TreatEvent(){
+ if(CheckGoodEvent()){
+ ReInit();
+ SetParticles();
+
+ TVector3 BeamDir = InitialConditions->GetBeamDirection();
+ IncidentTheta = BeamDir.Angle(ZDir);
+
+
+ std::vector<double> dE_IC(ZDD->GetICcounter(),0);
+
+ ///////////////////////////// Using Initial Conditions as CATS ///////////////
+ // This part seems to be working fine
+ double xinit = InitialConditions->GetIncidentPositionX();
+ double yinit = InitialConditions->GetIncidentPositionY();
+ double zinit = InitialConditions->GetIncidentPositionZ();
+ double ztarget = m_DetectorManager->GetTargetZ();
+ double coefficient = (ztarget - zinit)/BeamDir.Z(); // coeff such that k*zdir = (ztarget - zinit)
+ xtarget = xinit + coefficient*BeamDir.X();
+ ytarget = yinit + coefficient*BeamDir.y();
+ BeamImpact = TVector3(xtarget,ytarget,ztarget);
+ BeamDirection = TVector3(xtarget - xinit, ytarget - yinit, ztarget - zinit);
+
+ FinalBeamEnergy = Beam_Target.Slow(OriginalBeamEnergy,TargetThickness*0.5,BeamDirection.Angle(ZDir));
+ reaction->SetBeamEnergy(FinalBeamEnergy);
+
+
+ //////////////////////////// Treating DC Position ///////////////////////////
+
+ for(int i = 0; i < ZDD->DC_DetectorNumber.size(); i++){
+ if(ZDD->DC_DetectorNumber[i] == 1){
+ ZDD_DC_X = Drift_Speed*(ZDD->DC_DriftTime[i]*microsecond) - Drift_Chamber_Length/2;
+ }
+ else if(ZDD->DC_DetectorNumber[i] == 2){
+ ZDD_DC_Y = Drift_Speed*(ZDD->DC_DriftTime[i]*microsecond) - Drift_Chamber_Width/2;
+ }
+ }
+
+ // Faut encore réfléchir aux angles, y'a un truc pas net, l'angle après la cible est bcp trop large
+ TVector3 AfterTargetDir_X(ZDD_DC_X - xtarget, 0,ZDD_R + ZDD->Get_m_Drift_Chamber_Z()[0] + ZDD->Get_m_Drift_Chamber_Thickness()[0]*0.5 - ztarget);
+ ZDD_ThetaAfterTarget_X = AfterTargetDir_X.Angle(ZDir);
+ if(ZDD_DC_X < xtarget){
+ ZDD_ThetaAfterTarget_X *= -1;
+ }
+ TVector3 AfterTargetDir((ZDD_R + ZDD->Get_m_Drift_Chamber_Z()[1] + ZDD->Get_m_Drift_Chamber_Thickness()[1]*0.5 -ztarget)*tan(ZDD_ThetaAfterTarget_X) - xtarget,
+ ZDD_DC_Y - ytarget,
+ ZDD_R + ZDD->Get_m_Drift_Chamber_Z()[1] + ZDD->Get_m_Drift_Chamber_Thickness()[1]*0.5 - ztarget);
+ ZDD_ThetaAfterTarget = AfterTargetDir.Angle(ZDir);
+ /*if(tan(ZDD_ThetaAfterTarget_X) < 0){
+ ZDD_ThetaAfterTarget *= -1;
+ }*/
+ // Faire gaffe: se poser la question qui est X et qui est Y;
+
+ /////////////////////// Treating Energy ///////////////////////////////
+ ZDD_dE_tot = 0;
+ ZDD_E_tot = 0;
+ ZDD_E_Plastic = 0;
+ for(int i = 0; i < ZDD->Plastic_DetectorNumber.size(); i++){
+ ZDD_E_Plastic+= ZDD->Plastic_Energy[i];
+ }
+
+ for(int i = 0; i < ZDD->IC_DetectorNumber.size(); i++){
+ dE_IC[ZDD->IC_DetectorNumber[i]]+= ZDD->IC_Energy[i];
+ }
+
+ ////////////////////////// Note on how to compute full energy recovery: Plastic, then 1 plan of Kapton, then gas, then Mylar,
+ ////////////////////////// then all IC + Mylar, then again gas and Kapton, and Mylar at the entry and exit of each DC
+
+ //ZDD_ThetaAfterTarget = 0;
+ ZDD_E_tot = ZDD_E_Plastic;
+ ZDD_E_tot = Heavy_IC_Windows[1].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_Entry_Exit_Thickness()[1], ZDD_ThetaIC + ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_IC_Gas[1].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_m_Gas_Gap_Thickness()[1], ZDD_ThetaIC + ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_IC_Mylar[ZDD->Get_AC_Z().size()-1].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_AC_Thickness()[ZDD->Get_AC_Z().size()-1], ZDD_ThetaIC + ZDD_ThetaAfterTarget);
+ // Ok ça c'est bon
+ //
+ for(int i = ZDD->IC_Energy.size() -1; i >= 0; i--){
+ //if(ZDD->IC_Energy[i] > 5){
+ ZDD_E_tot+= ZDD->IC_Energy[i];
+ ZDD_dE_tot+= ZDD->IC_Energy[i];
+ ZDD_E_tot = Heavy_IC_Mylar[i+Nb_Mylar_Before_IC].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_AC_Thickness()[i+Nb_Mylar_Before_IC], ZDD_ThetaIC + ZDD_ThetaAfterTarget);
+ //}
+ }
+ ZDD_E_tot = Heavy_IC_Gas[0].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_m_Gas_Gap_Thickness()[0], ZDD_ThetaIC + ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_IC_Windows[0].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_Entry_Exit_Thickness()[0], ZDD_ThetaIC + ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_IC_Mylar[3].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_AC_Thickness()[3], ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_DC_Gas[1].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_m_Drift_Chamber_Thickness()[1], ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_IC_Mylar[2].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_AC_Thickness()[2], ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_IC_Mylar[1].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_AC_Thickness()[1], ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_DC_Gas[0].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_m_Drift_Chamber_Thickness()[0], ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_IC_Mylar[0].EvaluateInitialEnergy(ZDD_E_tot, ZDD->Get_AC_Thickness()[0], ZDD_ThetaAfterTarget);
+ ZDD_E_tot = Heavy_Target.EvaluateInitialEnergy(ZDD_E_tot, TargetThickness/2 , ZDD_ThetaAfterTarget);
+
+ //ZDD_E_tot = Beam_Target.EvaluateInitialEnergy(ZDD_E_tot, TargetThickness/2 , ZDD_ThetaAfterTarget);
+
+ // We reproduce the succession of materials crossed
+
+ //////////////////// MUST2 Part ////////////////////
+
+ for(unsigned int countMust2 = 0 ; countMust2 < M2->Si_E.size() ; countMust2++){
+ /************************************************/
+ //Part 0 : Get the usefull Data
+ // MUST2
+ int TelescopeNumber = M2->TelescopeNumber[countMust2];
+
+ /************************************************/
+ // Part 1 : Impact Angle
+ ThetaM2Surface = 0;
+ ThetaNormalTarget = 0;
+ TVector3 HitDirection = M2 -> GetPositionOfInteraction(countMust2) - BeamImpact ;
+ M2_ThetaLab = HitDirection.Angle( BeamDirection );
+
+ M2_X = M2 -> GetPositionOfInteraction(countMust2).X();
+ M2_Y = M2 -> GetPositionOfInteraction(countMust2).Y();
+ M2_Z = M2 -> GetPositionOfInteraction(countMust2).Z();
+
+ ThetaM2Surface = HitDirection.Angle(- M2 -> GetTelescopeNormal(countMust2) );
+ ThetaNormalTarget = HitDirection.Angle( TVector3(0,0,1) ) ;
+
+ /************************************************/
+
+ /************************************************/
+ // Part 2 : Impact Energy
+ Energy = M2_ELab = 0;
+ Si_E_M2 = M2->Si_E[countMust2];
+ CsI_E_M2= M2->CsI_E[countMust2];
+
+ // if CsI
+ if(CsI_E_M2>0 ){
+ // The energy in CsI is calculate form dE/dx Table because
+ Energy = CsI_E_M2;
+ Energy = LightAl.EvaluateInitialEnergy( Energy ,0.4*micrometer , ThetaM2Surface);
+ Energy+=Si_E_M2;
+ }
+
+ else
+ Energy = Si_E_M2;
+
+
+ // Evaluate energy using the thickness
+ M2_ELab = LightAl.EvaluateInitialEnergy( Energy ,0.4*micrometer , ThetaM2Surface);
+ // Target Correction
+ M2_ELab = LightTarget.EvaluateInitialEnergy( M2_ELab ,TargetThickness*0.5, ThetaNormalTarget);
+ M2_dE = Si_E_M2;
+
+ /************************************************/
+
+ /************************************************/
+ // Part 3 : Excitation Energy Calculation
+ //std::cout << reaction->GetBeamEnergy() << std::endl;
+ M2_Ex = reaction->ReconstructRelativistic( M2_ELab , M2_ThetaLab );
+ reaction->SetBeamEnergy(InitialConditions->GetIncidentFinalKineticEnergy());
+ M2_ExNoBeam=reaction->ReconstructRelativistic( M2_ELab , M2_ThetaLab );
+ reaction->SetBeamEnergy(FinalBeamEnergy);
+ M2_ExNoProton = reaction->ReconstructRelativistic( ReactionConditions->GetKineticEnergy(0) , ReactionConditions->GetParticleDirection(0).Angle(TVector3(0,0,1)) );
+ M2_ThetaLab=M2_ThetaLab/deg;
+
+ /************************************************/
+
+ /************************************************/
+ // Part 4 : Theta CM Calculation
+ M2_ThetaCM = reaction->EnergyLabToThetaCM( M2_ELab , M2_ThetaLab)/deg;
+ /************************************************/
+ }//end loop MUST2
+
+
+
+ }
+ //ZDD->Clear();
+}
+
+void Analysis::SetParticles(){
+ BeamPart->SetUp(BeamName);
+ BeamPart->SetKineticEnergy(InitialConditions->GetIncidentInitialKineticEnergy());
+ Beta = BeamPart->GetBeta();
+ Gamma = BeamPart->GetGamma();
+ Velocity = BeamPart->GetVelocity();
+}
+
+bool Analysis::CheckGoodEvent(){
+
+ return //CheckIC() && CheckPlastics() && CheckDC();
+ true;
+}
+
+bool Analysis::CheckIC(){
+
+ std::vector<int> IC_vec(ZDD->GetICcounter(), 0);
+ bool result = true;
+ for(int i = 0; i < ZDD->IC_DetectorNumber.size(); i++){
+ if(ZDD->IC_Energy[i] > 0){
+ IC_vec[ZDD->IC_DetectorNumber[i]]++;
+ }
+
+ }
+ for(int i = 0; i < ZDD->GetICcounter(); i++){
+ if(IC_vec[i] != 1){
+ result = false;
+ }
+ }
+ return result;
+}
+bool Analysis::CheckPlastics(){
+ int counter = 0;
+ for(int i = 0; i < ZDD->Plastic_DetectorNumber.size(); i++){
+ if(ZDD->Plastic_Energy[i] > 0){
+ counter++;
+ }
+ }
+ if(counter == 1){
+ return true;
+ }
+ else{
+ return false;
+ }
+}
+
+bool Analysis::CheckDC(){
+ int DC1_counter = 0, DC2_counter = 0;
+ for(int i = 0; i < ZDD->DC_DetectorNumber.size(); i++){
+ if(ZDD->DC_DetectorNumber[i] == 1){
+ DC1_counter++;
+ }
+ else if(ZDD->DC_DetectorNumber[i] == 2){
+ DC2_counter++;
+ }
+ }
+ return DC2_counter == 1 && DC1_counter == 1;
+}
+
+void Analysis::InitOutputBranch() {
+ /////////////////////// ZDD related branch ////////////////////////////////
+ RootOutput::getInstance()->GetTree()->Branch("ZDD_E_tot",&ZDD_E_tot,"ZDD_E_tot/D");
+ //RootOutput::getInstance()->GetTree()->Branch("ZDD_IC_Energy",&ZDD_IC_Energy,"ZDD_IC_Energy[10]/D");
+ RootOutput::getInstance()->GetTree()->Branch("ZDD_dE_tot",&ZDD_dE_tot,"ZDD_dE_tot/D");
+ RootOutput::getInstance()->GetTree()->Branch("ZDD_E_Plastic",&ZDD_E_Plastic,"ZDD_E_Plastic/D");
+ RootOutput::getInstance()->GetTree()->Branch("ZDD_DC_X",&ZDD_DC_X,"ZDD_DC_X/D");
+ RootOutput::getInstance()->GetTree()->Branch("ZDD_DC_Y",&ZDD_DC_Y,"ZDD_DC_Y/D");
+ RootOutput::getInstance()->GetTree()->Branch("xtarget",&xtarget,"xtarget/D");
+ RootOutput::getInstance()->GetTree()->Branch("ytarget",&ytarget,"ytarget/D");
+ RootOutput::getInstance()->GetTree()->Branch("Beta",&Beta,"Beta/D");
+ RootOutput::getInstance()->GetTree()->Branch("ZDD_ThetaAfterTarget",&ZDD_ThetaAfterTarget,"ZDD_ThetaAfterTarget/D");
+ RootOutput::getInstance()->GetTree()->Branch("ZDD_ThetaAfterTarget_X",&ZDD_ThetaAfterTarget_X,"ZDD_ThetaAfterTarget_X/D");
+ RootOutput:: getInstance()->GetTree()->Branch("InitialConditions",&InitialConditions);
+
+ //////////////////////// MUST related branch ////////////////////////////////
+ RootOutput::getInstance()->GetTree()->Branch("M2_Ex",&M2_Ex,"M2_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("M2_ELab",&M2_ELab,"M2_ELab/D");
+ RootOutput::getInstance()->GetTree()->Branch("M2_ThetaLab",&M2_ThetaLab,"M2_ThetaLab/D");
+ RootOutput::getInstance()->GetTree()->Branch("M2_ThetaCM",&M2_ThetaCM,"M2_ThetaCM/D");
+ RootOutput::getInstance()->GetTree()->Branch("M2_X",&M2_X,"M2_X/D");
+ RootOutput::getInstance()->GetTree()->Branch("M2_Y",&M2_Y,"M2_Y/D");
+ RootOutput::getInstance()->GetTree()->Branch("M2_Z",&M2_Z,"M2_Z/D");
+ RootOutput::getInstance()->GetTree()->Branch("M2_dE",&M2_dE,"M2_dE/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::ReInit(){
+ ZDD_E_tot = -1000;
+ ZDD_E_Plastic = -1000;
+ ZDD_dE_tot= -1000;
+ ZDD_DC_X = -1000;
+ ZDD_DC_Y = -1000;
+ xtarget = -1000;
+ ytarget = -1000;
+ Beta = -1000;
+
+ M2_Ex = -1000 ;
+ //ExNoBeam=ExNoProton=-1000;
+ //EDC= -1000;
+ M2_ELab = -1000;
+ //BeamEnergy = -1000;
+ M2_ThetaLab = -1000;
+ M2_ThetaCM = -1000;
+ M2_X = -1000;
+ M2_Y = -1000;
+ M2_Z = -1000;
+ M2_dE= -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/MUST_AND_ZDD/Analysis.h b/Projects/MUST_AND_ZDD/Analysis.h
new file mode 100644
index 0000000000000000000000000000000000000000..2a40fd9c58cd5ca87afd00a7518a72869b495e13
--- /dev/null
+++ b/Projects/MUST_AND_ZDD/Analysis.h
@@ -0,0 +1,146 @@
+#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 MUST_AND_ZDD analysis project *
+ * *
+ *---------------------------------------------------------------------------*
+ * Comment: *
+ * *
+ *****************************************************************************/
+
+#include "TMust2Physics.h"
+#include"NPVAnalysis.h"
+#include"TZDDPhysics.h"
+#include"NPEnergyLoss.h"
+#include"NPFunction.h"
+#include"NPReaction.h"
+#include"NPOptionManager.h"
+#include"RootInput.h"
+#include"TInitialConditions.h"
+#include "TReactionConditions.h"
+#include"NPParticle.h"
+#include"NPBeam.h"
+class Analysis: public NPL::VAnalysis{
+ public:
+ Analysis();
+ ~Analysis();
+
+ public:
+ void SetParticles();
+ void SetEnergyLoss();
+ void Init();
+ void TreatEvent();
+ void End();
+ void ReInit();
+ void InitOutputBranch();
+ void InitInputBranch();
+ bool CheckGoodEvent();
+ // CheckIC is checking that all IC are crossed (E>0) not less and not more than 1 time
+ bool CheckIC();
+ // CheckPlastics is checking that only and at leat 1 plastic is hit
+ bool CheckPlastics();
+ // CHecking multiplicity 1 in DC
+ bool CheckDC();
+
+ static NPL::VAnalysis* Construct();
+
+ private:
+
+ /// Currently only treating multiplicity 1 events
+ // ZDD info
+ double ZDD_DC_X;
+ double ZDD_DC_Y;
+ double ZDD_ThetaIC = 30/deg;
+ double ZDD_ThetaAfterTarget;
+ double ZDD_ThetaAfterTarget_X;
+ double ZDD_E_tot;
+ double ZDD_E_Plastic;
+ double ZDD_dE_tot;
+
+ // MUST2 info
+ double M2_Ex;
+ double M2_ExNoBeam;
+ double M2_ExNoProton;
+ double M2_EDC;
+ double M2_ELab;
+ double M2_ThetaLab;
+ double M2_ThetaCM;
+ double M2_X;
+ double M2_Y;
+ double M2_Z;
+ double M2_dE;
+
+ double OriginalBeamEnergy ; // AMEV
+ double FinalBeamEnergy;
+
+
+ double xtarget;
+ double ytarget;
+ double IncidentTheta = 0;
+ int DetectorNumber ;
+ double ThetaNormalTarget;
+ double ThetaM2Surface ;
+ double ThetaMGSurface ;
+ double Si_E_M2 ;
+ double CsI_E_M2 ;
+ double Energy ;
+ double BeamEnergy;
+ double ThetaGDSurface ;
+
+ double Beta;
+ double Gamma;
+ double Velocity;
+
+ double Drift_Speed;
+ double TargetThickness;
+ double CorrectedBeamEnergy;
+ std::vector<double> IC_Energy;
+
+ TVector3 ZDir;
+ TVector3 BeamDirection;
+ TVector3 BeamImpact;
+
+
+
+ NPL::Reaction* reaction = new Reaction;
+ NPL::EnergyLoss Beam_Target;
+ NPL::EnergyLoss Heavy_Target;
+ NPL::EnergyLoss LightTarget;
+ std::vector<NPL::EnergyLoss> Heavy_IC_Gas;
+ std::vector<NPL::EnergyLoss> Heavy_IC_Windows;
+ std::vector<NPL::EnergyLoss> Heavy_IC_Mylar;
+ std::vector<NPL::EnergyLoss> Heavy_DC_Gas;
+ NPL::EnergyLoss LightAl;
+ NPL::EnergyLoss LightSi;
+
+ string BeamName = "48Cr";
+
+ int Nb_Mylar_After_IC = 1;
+ int Nb_Mylar_Before_IC = 4;
+ double Drift_Chamber_Length = 400*mm;
+ double Drift_Chamber_Width = 400*mm;
+ double ZDD_R = 1*m;
+
+ private:
+ TMust2Physics* M2;
+ TZDDPhysics* ZDD;
+ TInitialConditions* InitialConditions;
+ TReactionConditions* ReactionConditions;
+
+ Beam* BeamPart = new Beam;
+ Particle* HeavyEjectile = new Particle;
+};
+#endif
\ No newline at end of file
diff --git a/Projects/MUST_AND_ZDD/CMakeLists.txt b/Projects/MUST_AND_ZDD/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..22c74affdfc45019bdda2594f8439c52d4ab97ec
--- /dev/null
+++ b/Projects/MUST_AND_ZDD/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/MUST_AND_ZDD/MUST_AND_ZDD.detector b/Projects/MUST_AND_ZDD/MUST_AND_ZDD.detector
new file mode 100644
index 0000000000000000000000000000000000000000..6b20f9a778b2e9a8534ace6210199e348c747d50
--- /dev/null
+++ b/Projects/MUST_AND_ZDD/MUST_AND_ZDD.detector
@@ -0,0 +1,312 @@
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+Target
+ THICKNESS= 53.5 micrometer
+ ANGLE= 0 deg
+ RADIUS= 15 mm
+ MATERIAL= CH2
+ ANGLE= 0 deg
+ X= 0 mm
+ Y= 0 mm
+ Z= 0 mm
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ZDD
+ R = 1 m
+ Theta = 0 deg
+
+
+ZDD AC
+ Z= 50 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=0 deg
+
+ZDD DC
+ Z= 50.002 mm
+ Thickness= 119.998 mm
+ Gas= iC4H10
+ Pressure= 0.006 bar
+ Temperature= 295 kelvin
+
+ZDD AC
+ Z= 170 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=0 deg
+
+ZDD AC
+ Z= 420 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=0 deg
+
+ZDD DC
+ Z= 420.002 mm
+ Thickness= 119.998 mm
+ Gas= iC4H10
+ Pressure= 0.006 bar
+ Temperature= 295 kelvin
+
+ZDD AC
+ Z= 540 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=0 deg
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 1
+ZDD AC
+ Z= 1500 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1500.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 2
+ZDD AC
+ Z= 1540 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1540.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 3
+ZDD AC
+ Z= 1580.00231 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1580.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 4
+ZDD AC
+ Z= 1620 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1620.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 5
+ZDD AC
+ Z= 1660 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1660.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 6
+ZDD AC
+ Z= 1700 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1700.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 7
+ZDD AC
+ Z= 1740 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1740.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 8
+ZDD AC
+ Z= 1780 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1780.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 9
+ZDD AC
+ Z= 1820 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1820.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 10
+ZDD AC
+ Z= 1860 mm
+ Thickness= 2.000 um
+ Material = Mylar
+ Theta=30 deg
+
+ZDD IC
+ Z= 1860.00231 mm
+ Thickness= 34.63769 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 11
+ZDD AC
+ Z= 1900 mm
+ Thickness= 1.000 um
+ Material = Mylar
+ Theta=30 deg
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Entry_Exit
+ZDD EntryExit
+ Z= 1960 mm
+ Thickness= 10.00 um
+ Material = Kapton
+
+ZDD EntryExit
+ Z= 1490 mm
+ Thickness= 10.00 um
+ Material = Kapton
+
+ZDD GasGap
+ Z= 1490.01155 mm
+ Thickness= 8.64845 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+ZDD GasGap
+ Z= 1900.001155 mm
+ Thickness= 43.28975 mm
+ Gas= CF4
+ Pressure= 0.5 bar
+ Temperature= 295 kelvin
+
+
+
+ZDD Plastic
+ Material= BC400
+ Width= 1000 mm
+ Length= 100 mm
+ Thickness= 30 mm
+ Pos= 0 -200 3000 mm
+
+ZDD Plastic
+ Material= BC400
+ Width= 1000 mm
+ Length= 100 mm
+ Thickness= 30 mm
+ Pos= 0 -100 3000 mm
+
+ZDD Plastic
+ Material= BC400
+ Width= 1000 mm
+ Length= 100 mm
+ Thickness= 30 mm
+ Pos= 0 0 3000 mm
+
+ZDD Plastic
+ Material= BC400
+ Width= 1000 mm
+ Length= 100 mm
+ Thickness= 30 mm
+ Pos= 0 100 3000 mm
+
+ZDD Plastic
+ Material= BC400
+ Width= 1000 mm
+ Length= 100 mm
+ Thickness= 30 mm
+ Pos= 0 200 3000 mm
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%%%%%%% Telescope 1 %%%%%%%
+M2Telescope
+ X1_Y1= 10.851 105.032 298.604 mm
+ X1_Y128= 22.798 9.835 328.389 mm
+ X128_Y1= 104.089 105.032 261.204 mm
+ X128_Y128= 116.036 9.835 290.089 mm
+ SI= 1
+ SILI= 0
+ CSI= 1
+ VIS= all
+
+%%%%%%% Telescope 2 %%%%%%%
+M2Telescope
+ X1_Y1= -116.036 9.835 290.089 mm
+ X1_Y128= -22.798 9.835 328.389 mm
+ X128_Y1= -104.089 105.032 261.204 mm
+ X128_Y128= -10.851 105.032 298.604 mm
+ SI= 1
+ SILI= 0
+ CSI= 1
+ VIS= all
+
+
+%%%%%%% Telescope 3 %%%%%%%
+M2Telescope
+ X1_Y1= -10.851 -105.032 298.604 mm
+ X1_Y128= -22.798 -9.835 328.389 mm
+ X128_Y1= -104.089 -105.032 261.204 mm
+ X128_Y128= -116.036 -9.835 290.089 mm
+ SI= 1
+ SILI= 0
+ CSI= 1
+ VIS= all
+
+
+
+%%%%%%% Telescope 4 %%%%%%%
+M2Telescope
+ X1_Y1= 116.036 -9.835 290.089 mm
+ X1_Y128= 22.798 -9.835 328.389 mm
+ X128_Y1= 104.089 -105.032 261.204 mm
+ X128_Y128= 10.851 -105.032 298.604 mm
+ SI= 1
+ SILI= 0
+ CSI= 1
+ VIS= all
diff --git a/Projects/e793s/Analysis.cxx b/Projects/e793s/Analysis.cxx
index c04cbfdd958e03e3fbb5f0552925b3f42c6cb1d5..aa1eabc63a3665dd96697432752483cd743fcdbb 100755
--- a/Projects/e793s/Analysis.cxx
+++ b/Projects/e793s/Analysis.cxx
@@ -227,9 +227,11 @@ void Analysis::TreatEvent(){
/**/
if(isSim && !isPhaseSpace){
- if(M2->Si_E[countMust2]>0 && M2->CsI_E[countMust2]<=0){
- /* IS a count in DSSD, but NOT in CsI to exclude punch through */
- /* May need further gating in Si_E & _T */
+ if(M2->TelescopeNumber[countMust2]<5){
+ if(M2->Si_E[countMust2]>0 && // DSSD count
+ M2->CsI_E[countMust2]<=0 && // No CsI count
+ M2->Si_T[countMust2]<460 // Triton kinematic line, not punch through
+ ){
ThetaCM_detected_MM->Fill(ReactionConditions->GetThetaCM());
ThetaLab_detected_MM->Fill(ReactionConditions->GetTheta(0));
@@ -237,6 +239,15 @@ void Analysis::TreatEvent(){
ThetaCM_detected_MMX[MMX]->Fill(ReactionConditions->GetThetaCM());
ThetaLab_detected_MMX[MMX]->Fill(ReactionConditions->GetTheta(0));
}
+ } else {
+ //No triton requirement for MM5
+ ThetaCM_detected_MM->Fill(ReactionConditions->GetThetaCM());
+ ThetaLab_detected_MM->Fill(ReactionConditions->GetTheta(0));
+
+ int MMX = TelescopeNumber-1;
+ ThetaCM_detected_MMX[MMX]->Fill(ReactionConditions->GetThetaCM());
+ ThetaLab_detected_MMX[MMX]->Fill(ReactionConditions->GetTheta(0));
+ }
}
/**/
diff --git a/Projects/e793s/Analysis.h b/Projects/e793s/Analysis.h
index 071c320aaa57eb6023ed4aa93f71ad0d8533c6e8..150d05057f51ff658658845895a4bf71bec4fe36 100755
--- a/Projects/e793s/Analysis.h
+++ b/Projects/e793s/Analysis.h
@@ -42,7 +42,7 @@
#include <TMath.h>
#include <bitset>
-int NumThetaAngleBins = 900;// 180 = 1 deg, 360 = 0.5 deg, 900 = 0.2 deg, 1800 = 0.1 deg
+int NumThetaAngleBins = 1800;// 180 = 1 deg, 360 = 0.5 deg, 900 = 0.2 deg, 1800 = 0.1 deg
auto ThetaCM_emmitted = new TH1F("ThetaCM_emmitted","ThetaCM_emmitted",NumThetaAngleBins,0,180);
auto ThetaCM_detected_MM = new TH1F("ThetaCM_detected_MM","ThetaCM_detected_MM",NumThetaAngleBins,0,180);
diff --git a/Projects/e793s/macro/CS2_dt.h b/Projects/e793s/macro/CS2_dt.h
index a479ea275bdea7a908d86d33d8c4736453ca1ba7..87614b9b695c87058152dcfd1b1090e8c67736ac 100644
--- a/Projects/e793s/macro/CS2_dt.h
+++ b/Projects/e793s/macro/CS2_dt.h
@@ -94,8 +94,9 @@ void CS_Diagnosis(){
void CS(){
/* Overload function */
cout << "- CS(stateE, stateSp, orb_l, orb_j, nodes) "<< endl;
- cout << "---- 1.945, p1/2 = CS(1.945, 1, 0, 0.5) "<< endl;
+ cout << "---- 1.945, s1/2 = CS(1.945, 1, 0, 0.5) "<< endl;
cout << "---- 1.945, d3/2 = CS(1.945, 1, 2, 1.5) "<< endl;
+ cout << "---- 0.691, f7/2 = CS(0.691, 4, 3, 3.5) "<< endl;
}
////////////////////////////////////////////////////////////////////////////////
@@ -536,23 +537,19 @@ vector<vector<double>> GetExpDiffCross(double Energy){
/* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! */
/* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! */
- // TEMPORARY FIX FOR EX:THETACM PROBLEM!!
- // Assuming angular bins are 2-20 in 1 degree sections...
- if(means_dt[indexE] > 4.0){
- numAngleBins-= 4; //Up to 16 deg
- }
- else if(means_dt[indexE] > 3.0){
- numAngleBins-= 5; //Up to 15 deg
- }
- else if(means_dt[indexE] > 2.0){
- numAngleBins-= 6; //Up to 14 deg
- }
- else if(means_dt[indexE] > 1.0){
- numAngleBins-= 8; //Up to 12 deg
- }
- else if(means_dt[indexE] > 0.0){
- numAngleBins-=12; //Up to 8 deg
- }
+ /* TEMPORARY FIX FOR EX:THETACM PROBLEM!! */
+ double thetaMax = 0.23*pow(means_dt[indexE],2)
+ + 0.76*means_dt[indexE]
+ + 9.01;
+ thetaMax = floor(thetaMax);
+
+ double thetaMin = 0.10*pow(means_dt[indexE],2)
+ + 0.06*means_dt[indexE]
+ + 1.58;
+ thetaMin = ceil(thetaMin);
+
+ numAngleBins = (int)thetaMax - (int)thetaMin;
+ firstAngle = (int)thetaMin;
/* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! */
/* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! */
diff --git a/Projects/e793s/macro/DrawPlots.h b/Projects/e793s/macro/DrawPlots.h
index 4456c7e58f5d4cd6ac88822052ea0ede20f41593..44eac574080ab87a6e67a30d3edd047ada187414 100755
--- a/Projects/e793s/macro/DrawPlots.h
+++ b/Projects/e793s/macro/DrawPlots.h
@@ -1199,6 +1199,39 @@ void ExThetaLab(){
Sn->Draw();
}
+void ExThetaCM(){
+ string gate = timegate
+ + " && " + det_gate
+ + " && Ex@.size()==1";
+ if(reactionName=="47K(d,t)"){
+ gate = gate + " && cutTritons && cutTime";
+ }
+
+ TCanvas *diagnoseTheta = new TCanvas("diagnoseTheta","diagnoseTheta",1000,1000);
+ chain->Draw(
+ "Ex:ThetaCM>>thetaHistCM(360,0,180,180,-1,8)",
+ gate.c_str(), "colz");
+ TH1F* thetaHistCM = (TH1F*) gDirectory->Get("thetaHistCM");
+ thetaHistCM->GetXaxis()->SetTitle("#theta_{CM} [deg]");
+ thetaHistCM->GetYaxis()->SetTitle("Ex [MeV]");
+ thetaHistCM->SetTitle("Theta dependance testing");
+
+ diagnoseTheta->Update();
+
+ TLine *l0000 = new TLine(100., 0.000, 160., 0.000);
+ l0000->SetLineStyle(kDashed);
+ l0000->SetLineColor(kRed);
+ l0000->Draw();
+
+ if(reactionName=="47K(d,p)"){
+ TLine *Sn = new TLine(100., 4.644, 160., 4.644);
+ Sn->SetLineStyle(kDashed);
+ Sn->SetLineColor(kRed);
+ Sn->Draw();
+ }
+}
+
+
void ExThetaLab(double gamma, double width){
TCanvas *diagnoseTheta = new TCanvas("diagnoseTheta","diagnoseTheta",1000,1000);
string gate = timegate
@@ -2864,3 +2897,120 @@ void Figure_TopGamma_BottomParticle(double gammaBinWidth, double particleBinWidt
cFig_GGSP->Update();
}
+
+void Figure_SolidAngle_dt(){
+ //string histname = "SolidAngle_CM_MM";
+ string histname = "SolidAngle_Lab_MM";
+
+ TCanvas* canv = new TCanvas("canv","canv",1000,1000);
+ gStyle->SetPadLeftMargin(0.10);
+ gStyle->SetPadRightMargin(0.03);
+ gStyle->SetOptStat(0);
+
+ TFile *file6 = new TFile("./SolidAngle_HistFiles/SAHF_18Oct22_47Kdt_6000.root","READ");
+ TFile *file5 = new TFile("./SolidAngle_HistFiles/SAHF_18Oct22_47Kdt_5000.root","READ");
+ TFile *file4 = new TFile("./SolidAngle_HistFiles/SAHF_18Oct22_47Kdt_4000.root","READ");
+ TFile *file3 = new TFile("./SolidAngle_HistFiles/SAHF_18Oct22_47Kdt_3000.root","READ");
+ TFile *file2 = new TFile("./SolidAngle_HistFiles/SAHF_18Oct22_47Kdt_2000.root","READ");
+ TFile *file1 = new TFile("./SolidAngle_HistFiles/SAHF_18Oct22_47Kdt_1000.root","READ");
+ TFile *file0 = new TFile("./SolidAngle_HistFiles/SAHF_18Oct22_47Kdt_0000.root","READ");
+
+ TH1F* h6 = (TH1F*)file6->Get(histname.c_str());
+ //h6->GetXaxis()->SetTitle("#theta_{CM}");
+ h6->GetXaxis()->SetTitle("#theta_{Lab}");
+ h6->GetYaxis()->SetTitle("Solid Angle [(sr)?]");
+ h6->SetTitle("Simulated triton solid angle variation with 46K state energy");
+ h6->GetXaxis()->SetRangeUser(0.,30.);
+ h6->GetYaxis()->SetRangeUser(0.,0.002);
+ h6->SetLineColor(kCyan+4);
+ h6->SetFillColor(kCyan+4);
+ h6->SetLineWidth(2);
+ h6->SetFillStyle(3002);
+ h6->Draw("hist");
+
+ TH1F* h5 = (TH1F*)file5->Get(histname.c_str());
+ h5->SetLineColor(kCyan-1);
+ h5->SetFillColor(kCyan-1);
+ h5->SetLineWidth(2);
+ h5->SetFillStyle(3002);
+ h5->Draw("hist same");
+
+ TH1F* h4 = (TH1F*)file4->Get(histname.c_str());
+ h4->SetLineColor(kCyan-5);
+ h4->SetFillColor(kCyan-5);
+ h4->SetLineWidth(2);
+ h4->SetFillStyle(3002);
+ h4->Draw("hist same");
+
+ TH1F* h3 = (TH1F*)file3->Get(histname.c_str());
+ h3->SetLineColor(kCyan-8);
+ h3->SetFillColor(kCyan-8);
+ h3->SetLineWidth(2);
+ h3->SetFillStyle(3002);
+ h3->Draw("hist same");
+
+ TH1F* h2 = (TH1F*)file2->Get(histname.c_str());
+ h2->SetLineColor(kCyan-9);
+ h2->SetFillColor(kCyan-9);
+ h2->SetLineWidth(2);
+ h2->SetFillStyle(3002);
+ h2->Draw("hist same");
+
+ TH1F* h1 = (TH1F*)file1->Get(histname.c_str());
+ h1->SetLineColor(kCyan-7);
+ h1->SetFillColor(kCyan-7);
+ h1->SetLineWidth(2);
+ h1->SetFillStyle(3002);
+ h1->Draw("hist same");
+
+ TH1F* h0 = (TH1F*)file0->Get(histname.c_str());
+ h0->SetLineColor(kCyan-0);
+ h0->SetFillColor(kCyan-0);
+ h0->SetLineWidth(2);
+ h0->SetFillStyle(3002);
+ h0->Draw("hist same");
+
+ TLatex *n0 = new TLatex(.5,.5,"0.0 MeV");
+ n0->SetTextColor(kCyan-0);
+ n0->SetTextSize(0.05);
+ n0->SetX(22.);
+ n0->SetY(0.0018);
+ n0->Draw("same");
+ TLatex *n1 = new TLatex(.5,.5,"1.0 MeV");
+ n1->SetTextColor(kCyan-7);
+ n1->SetTextSize(0.05);
+ n1->SetX(22.);
+ n1->SetY(0.0017);
+ n1->Draw("same");
+ TLatex *n2 = new TLatex(.5,.5,"2.0 MeV");
+ n2->SetTextColor(kCyan-9);
+ n2->SetTextSize(0.05);
+ n2->SetX(22.);
+ n2->SetY(0.0016);
+ n2->Draw("same");
+ TLatex *n3 = new TLatex(.5,.5,"3.0 MeV");
+ n3->SetTextColor(kCyan-8);
+ n3->SetTextSize(0.05);
+ n3->SetX(22.);
+ n3->SetY(0.0015);
+ n3->Draw("same");
+ TLatex *n4 = new TLatex(.5,.5,"4.0 MeV");
+ n4->SetTextColor(kCyan-5);
+ n4->SetTextSize(0.05);
+ n4->SetX(22.);
+ n4->SetY(0.0014);
+ n4->Draw("same");
+ TLatex *n5 = new TLatex(.5,.5,"5.0 MeV");
+ n5->SetTextColor(kCyan-1);
+ n5->SetTextSize(0.05);
+ n5->SetX(22.);
+ n5->SetY(0.0013);
+ n5->Draw("same");
+ TLatex *n6 = new TLatex(.5,.5,"6.0 MeV");
+ n6->SetTextColor(kCyan+4);
+ n6->SetTextSize(0.05);
+ n6->SetX(22.);
+ n6->SetY(0.0012);
+ n6->Draw("same");
+
+}
diff --git a/Projects/e793s/macro/KnownPeakFitter.h b/Projects/e793s/macro/KnownPeakFitter.h
index 42b6a3ccd6cfe4d6ac0f5476bbf7cc417a2fba2d..f4f66e89c990396dd31411fc02525c27dca15d35 100644
--- a/Projects/e793s/macro/KnownPeakFitter.h
+++ b/Projects/e793s/macro/KnownPeakFitter.h
@@ -24,10 +24,17 @@ array<double,numPeaks> means = { 0.000,
5.82
};
-const int numPeaks_dt = 3;
+const int numPeaks_dt = 8;//9;
array<double,numPeaks_dt> means_dt = { 0.000,
+ //0.587, //from lit
+ 0.691, //from lit
+ //0.886, //from lit
+ 1.738, //from lit
1.945,
+ 2.26 , //from lit
+ 2.76 , //from lit
3.344,
+ 4.2
};
@@ -371,3 +378,9 @@ void FitKnownPeaks(TH1F* hist){
//Shell function to call Rtrn_Arry without writing vector<vector<double>> to screen
vector<vector<double>> shell = FitKnownPeaks_RtrnArry(hist,0.0);
}
+
+void FitKnownPeaks_dt(TH1F* hist){
+ //Shell function to call Rtrn_Arry without writing vector<vector<double>> to screen
+ vector<vector<double>> shell = FitKnownPeaks_dt_RtrnArry(hist,0.0);
+}
+