diff --git a/Examples/Example1/ShowResults.C b/Examples/Example1/ShowResults.C
index 3a279292678101949dcbcc2813d5f7aafac7fde0..eabef9bf3da7fce25623e7d096d46b3c72f27876 100644
--- a/Examples/Example1/ShowResults.C
+++ b/Examples/Example1/ShowResults.C
@@ -36,7 +36,7 @@ void ShowResults(){
// E-TOF
c1->cd(2);
- chain->Draw("-MUST2.Si_T:SSSD.Energy+MUST2.Si_E>>hIDT(1000,0,35,1000,-30,0)","MUST2.CsI_E<0 && MUST2.TelescopeNumber<5","colz");
+ chain->Draw("MUST2.Si_T:SSSD.Energy+MUST2.Si_E>>hIDT(1000,0,35,500,450,500)","MUST2.CsI_E<0 && MUST2.TelescopeNumber<5","colz");
ETOF->Draw("same");
diff --git a/Examples/Example1/configs/ConfigMust2.dat b/Examples/Example1/configs/ConfigMust2.dat
index bcbc385cb42bc8a1e67a654e1877530ce21787f9..bebb0bf9ea475fdc1d492417228af5e96dbe37c5 100755
--- a/Examples/Example1/configs/ConfigMust2.dat
+++ b/Examples/Example1/configs/ConfigMust2.dat
@@ -1,5 +1,5 @@
ConfigMust2
- MAX_STRIP_MULTIPLICITY 1
+ MAX_STRIP_MULTIPLICITY 100
STRIP_ENERGY_MATCHING_NUMBER_OF_SIGMA 5
STRIP_ENERGY_MATCHING_SIGMA 0.02
DISABLE_CHANNEL MM1STRY12
@@ -14,6 +14,4 @@ ConfigMust2
DISABLE_ALL MM6
DISABLE_ALL MM7
DISABLE_ALL MM8
- SI_X_E_RAW_THRESHOLD 0
- CSI_E_RAW_THRESHOLD 0
CSI_SIZE 256
diff --git a/Examples/Example1/cuts/ETOF.root b/Examples/Example1/cuts/ETOF.root
index 067d548e7c6f2d4adc224047d2a23700cc380032..a57897ebf177e383f164c5d592085d98cb81da22 100644
Binary files a/Examples/Example1/cuts/ETOF.root and b/Examples/Example1/cuts/ETOF.root differ
diff --git a/NPLib/Core/NPCore.cxx b/NPLib/Core/NPCore.cxx
index 04fda36c45f4739c38638f7ccaea90110643bef9..35c76377e65b37a4481fb21bf270a6f068d051c9 100644
--- a/NPLib/Core/NPCore.cxx
+++ b/NPLib/Core/NPCore.cxx
@@ -65,3 +65,25 @@ namespace NPL{
std::string NPL::itoa(const int& i){
return NPL::itoa_array[i];
}
+////////////////////////////////////////////////////////////////////////////////
+unsigned int NPL::EnergyToADC(const double& E, const double& Emin, const double& Emax, const int& Pedestal, const int& ADCMax){
+ double Crange = ADCMax - Pedestal;
+ double Erange = Emax - Emin;
+
+ // Standard case where larger ADC channel means larger Energy
+ if(Crange > 0){
+ double answer = Crange*E/Erange + Pedestal;
+ if(answer > ADCMax)
+ answer = ADCMax;
+ return answer;
+ }
+ // Reverse coding case where larger ADC channel means smaller Energy
+ else{
+ double answer = Crange*E/Erange + Pedestal;
+ if(answer < ADCMax)
+ answer = ADCMax;
+ return answer;
+ }
+
+ return -1000;
+}
diff --git a/NPLib/Core/NPCore.h b/NPLib/Core/NPCore.h
index 9f14ffadd2b65bd0d096f8b1241805930713dc1e..0cf9ac1953b8ae64e25c10de220de77a1b24f200 100644
--- a/NPLib/Core/NPCore.h
+++ b/NPLib/Core/NPCore.h
@@ -28,5 +28,7 @@ namespace NPL{
void SendWarning(std::string Class, std::string Warning);
void SendInformation(std::string Class, std::string Information);
void SendErrorAndExit(std::string Class,std::string Error);
+ // For use to simulate electronic
+ unsigned int EnergyToADC(const double& E, const double& Emin, const double& Emax, const int& Pedestal, const int& ADCMax);
}
#endif
diff --git a/NPLib/Detectors/MUST2/TMust2Physics.cxx b/NPLib/Detectors/MUST2/TMust2Physics.cxx
index 7c5f7d735f685f378ddb4ca5662658a6d62a9f26..41a672f0f3d4479f242dd5079d19a97ad91e4982 100644
--- a/NPLib/Detectors/MUST2/TMust2Physics.cxx
+++ b/NPLib/Detectors/MUST2/TMust2Physics.cxx
@@ -896,24 +896,32 @@ map< string , TH1*> TMust2Physics::GetSpectra() {
void TMust2Physics::AddParameterToCalibrationManager()
{
CalibrationManager* Cal = CalibrationManager::getInstance();
+ // Good for simulation, close to typical values
+ vector<double> standardX = {-63, 63./8192.};
+ vector<double> standardY = {63, -63./8192.};
+ vector<double> standardCsI = {-250, 250./8192.};
+ vector<double> standardSiLi = {-250, 250./8192.};
+ vector<double> standardT = {-1000, 1000./8192.};
+
+
for(int i = 0 ; i < m_NumberOfTelescope ; ++i){
for( int j = 0 ; j < 128 ; ++j){
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_E") ;
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_E") ;
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_T") ;
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_T") ;
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_E",standardX) ;
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_E",standardY) ;
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_Si_X"+NPL::itoa(j+1)+"_T",standardT ) ;
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_Si_Y"+NPL::itoa(j+1)+"_T",standardT ) ;
}
for( int j = 0 ; j < 16 ; ++j){
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_E") ;
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_T") ;
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_E",standardSiLi) ;
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_SiLi"+NPL::itoa(j+1)+"_T",standardT ) ;
}
for( int j = 0 ; j < 16 ; ++j){
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_E") ;
- Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_T") ;
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_E","MUST2_T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_E",standardCsI);
+ Cal->AddParameter("MUST2", "T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_T","MUST2_T"+NPL::itoa(i+1)+"_CsI"+NPL::itoa(j+1)+"_T",standardT ) ;
}
}
diff --git a/NPSimulation/Detectors/MUST2/MUST2Array.cc b/NPSimulation/Detectors/MUST2/MUST2Array.cc
index 6796a6dfac3206508ac163a23905c90f850df6df..950046182d705c9fb0bafc476bf1046bda5d1340 100644
--- a/NPSimulation/Detectors/MUST2/MUST2Array.cc
+++ b/NPSimulation/Detectors/MUST2/MUST2Array.cc
@@ -44,7 +44,8 @@
#include "SiliconScorers.hh"
#include "CalorimeterScorers.hh"
#include "NPOptionManager.h"
-
+// NPL
+#include "NPCore.h"
//ROOT
#include "RootOutput.h"
@@ -758,9 +759,10 @@ void MUST2Array::ReadSensitive(const G4Event* event){
if(r==0){
trig.insert(detectorNbr);
// Energy
- m_Event->SetStripXE(detectorNbr,b+1,energyX);
+ m_Event->SetStripXE(detectorNbr,b+1,NPL::EnergyToADC(energyX,0,63,8192,16384));
// Time
- m_Event->SetStripXT(detectorNbr,b+1,RandGauss::shoot(time, ResoTimeMust));
+ double timeX = TimeOffset - RandGauss::shoot(time, ResoTimeMust);
+ m_Event->SetStripXT(detectorNbr,b+1,NPL::EnergyToADC(timeX,0,1000,8192,16384));
}
else{ // Interstrip X, keep maximum shared energy
double rand = G4UniformRand();
@@ -769,9 +771,10 @@ void MUST2Array::ReadSensitive(const G4Event* event){
if(energyX>0.1*keV){
trig.insert(detectorNbr);
// Energy
- m_Event->SetStripXE(detectorNbr,b+1,energyX) ;
+ m_Event->SetStripXE(detectorNbr,b+1,NPL::EnergyToADC(energyX,0,63,8192,16384)) ;
// Time
- m_Event->SetStripXT(detectorNbr,b+1,RandGauss::shoot(time, ResoTimeMust));
+ double timeX = TimeOffset - RandGauss::shoot(time, ResoTimeMust);
+ m_Event->SetStripXT(detectorNbr,b+1,NPL::EnergyToADC(timeX,0,1000,8192,16384));
}
}
else{
@@ -780,10 +783,11 @@ void MUST2Array::ReadSensitive(const G4Event* event){
trig.insert(detectorNbr);
// Energy
- m_Event->SetStripXE(detectorNbr,g+1,energyX) ;
+ m_Event->SetStripXE(detectorNbr,g+1,NPL::EnergyToADC(energyX,0,63,8192,16384)) ;
// Time
- m_Event->SetStripXT(detectorNbr,g+1,RandGauss::shoot(time, ResoTimeMust));
- }
+ double timeX = TimeOffset - RandGauss::shoot(time, ResoTimeMust);
+ m_Event->SetStripXT(detectorNbr,b+1,NPL::EnergyToADC(timeX,0,1000,8192,16384));
+ }
}
}
}
@@ -799,10 +803,11 @@ void MUST2Array::ReadSensitive(const G4Event* event){
if(r==0){
trig.insert(detectorNbr);
// Energy
- m_Event->SetStripYE(detectorNbr,b+1,energyY);
+ m_Event->SetStripYE(detectorNbr,b+1,NPL::EnergyToADC(energyY,0,63,8192,0));
// Time
- m_Event->SetStripYT(detectorNbr,b+1,RandGauss::shoot(time, ResoTimeMust));
- }
+ double timeY = TimeOffset - RandGauss::shoot(time, ResoTimeMust);
+ m_Event->SetStripYT(detectorNbr,b+1,NPL::EnergyToADC(timeY,0,1000,8192,16384));
+ }
else{ // Interstrip Y, keep both strip with shared energy
double rand = G4UniformRand();
double energyY1 = rand*energyY;
@@ -810,19 +815,22 @@ void MUST2Array::ReadSensitive(const G4Event* event){
trig.insert(detectorNbr);
// Energy
- m_Event->SetStripYE(detectorNbr,b+1,energyY1);
+ m_Event->SetStripYE(detectorNbr,b+1,NPL::EnergyToADC(energyY1,0,63,8192,0));
+
// Time
- m_Event->SetStripYT(detectorNbr,b+1,RandGauss::shoot(time, ResoTimeMust));
- }
+ double timeY = TimeOffset - RandGauss::shoot(time, ResoTimeMust);
+ m_Event->SetStripYT(detectorNbr,b+1,NPL::EnergyToADC(timeY,0,1000,8192,16384));
+ }
if(energyY1>0.1*keV){
trig.insert(detectorNbr);
double energyY2 = (1-rand)*energyY;
// Energy
- m_Event->SetStripYE(detectorNbr,g+1,energyY2);
+ m_Event->SetStripYE(detectorNbr,g+1,NPL::EnergyToADC(energyY2,0,63,8192,0));
// Time
- m_Event->SetStripYT(detectorNbr,g+1,RandGauss::shoot(time, ResoTimeMust));
- }
+ double timeY = TimeOffset - RandGauss::shoot(time, ResoTimeMust);
+ m_Event->SetStripYT(detectorNbr,b+1,NPL::EnergyToADC(timeY,0,1000,8192,16384));
+ }
}
}
@@ -859,8 +867,10 @@ void MUST2Array::ReadSensitive(const G4Event* event){
for(SiLi_itr = SiLiHitMap->GetMap()->begin(); SiLi_itr!=SiLiHitMap->GetMap()->end() ; SiLi_itr++){
G4double* Info = *(SiLi_itr->second);
if(Info[7]==*itr){//matching telescope number
- m_Event->SetSiLiE(Info[7],Info[8],RandGauss::shoot(Info[0],ResoSiLi));
- m_Event->SetSiLiT(Info[7],Info[8],RandGauss::shoot(Info[1],ResoTimeMust));
+ double ESiLi = RandGauss::shoot(Info[0],ResoSiLi);
+ m_Event->SetSiLiE(Info[7],Info[8],NPL::EnergyToADC(ESiLi,0,250,8192,16384));
+ double timeSiLi = RandGauss::shoot(Info[1],ResoTimeMust);
+ m_Event->SetSiLiT(Info[7],Info[8],NPL::EnergyToADC(timeSiLi,0,1000,16384,8192));
}
}
}
@@ -871,8 +881,10 @@ void MUST2Array::ReadSensitive(const G4Event* event){
G4double* Info = *(CsI_itr->second);
if(Info[7]==*itr){//matching telescope number
- m_Event->SetCsIE(Info[7],Info[8],RandGauss::shoot(Info[0],ResoCsI));
- m_Event->SetCsIT(Info[7],Info[8],RandGauss::shoot(Info[1],ResoTimeMust));
+ double ECsI = RandGauss::shoot(Info[0],ResoCsI);
+ m_Event->SetCsIE(Info[7],Info[8],NPL::EnergyToADC(ECsI,0,250,8192,16384));
+ double timeCsI = RandGauss::shoot(Info[1],ResoTimeMust);
+ m_Event->SetCsIT(Info[7],Info[8],NPL::EnergyToADC(timeCsI,0,1000,16384,8192));
}
}
}
diff --git a/NPSimulation/Detectors/MUST2/MUST2Array.hh b/NPSimulation/Detectors/MUST2/MUST2Array.hh
index ad42412f731a20931131f46bbdb6638afcec0899..883b62a014f0a5f593fcf54d3be2c750052f6e86 100644
--- a/NPSimulation/Detectors/MUST2/MUST2Array.hh
+++ b/NPSimulation/Detectors/MUST2/MUST2Array.hh
@@ -39,7 +39,7 @@ namespace MUST2
const G4double ResoSiLi = 0.028 ;// = 130 keV of resolution // Unit is MeV/2.35
const G4double ResoCsI = 0.08 ;// = 188 kev of resolution // Unit is MeV/2.35
const G4double ResoStrip = 0.0149 ;// 0.0223 = 52keV of Resolution // Unit is MeV/2.35 14.861996
-
+ const G4double TimeOffset = 500 ;// 500 ns stop
// Geometry
const G4double FaceFront = 11.*cm ;
const G4double FaceBack = 16.5*cm ;