Commit 11e50718 authored by Elidiano Tronchin's avatar Elidiano Tronchin
Browse files

Merge branch 'NPTool.2.dev' of https://github.com/adrien-matta/nptool into NPTool.2.dev

parents d199ec9e d95a3e0b
add_custom_command(OUTPUT TExlPhysicsDict.cxx COMMAND ../../scripts/build_dict.sh TExlPhysics.h TExlPhysicsDict.cxx TExlPhysics.rootmap libNPEXL.dylib DEPENDS TExlPhysics.h)
add_custom_command(OUTPUT TExlDataDict.cxx COMMAND ../../scripts/build_dict.sh TExlData.h TExlDataDict.cxx TExlData.rootmap libNPEXL.dylib DEPENDS TExlData.h)
add_library(NPEXL SHARED TExlData.cxx TExlPhysics.cxx TExlDataDict.cxx TExlPhysicsDict.cxx )
add_custom_command(OUTPUT TExlSpectraDict.cxx COMMAND ../../scripts/build_dict.sh TExlSpectra.h TExlSpectraDict.cxx TExlSpectra.rootmap libNPEXL.dylib DEPENDS TExlSpectra.h)
add_library(NPEXL SHARED TExlData.cxx TExlPhysics.cxx TExlSpectra.cxx TExlDataDict.cxx TExlPhysicsDict.cxx TExlSpectraDict.cxx)
target_link_libraries(NPEXL ${ROOT_LIBRARIES} NPCore)
install(FILES TExlData.h TExlPhysics.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY})
install(FILES TExlData.h TExlPhysics.h TExlSpectra.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY})
......@@ -22,39 +22,71 @@
#include <iostream>
#include "TExlData.h"
ClassImp(TExlData)
TExlData::TExlData()
{
}
TExlData::~TExlData()
{
}
void TExlData::Clear()
{
fExl_Energy.clear();
fExl_Time.clear();
fExl_Number.clear();
fExl_Crystal.clear();
}
fExlE_Energy.clear();
fExlE_Number.clear();
fExlE_Crystal.clear();
fExlT_Time.clear();
fExlT_Number.clear();
fExlT_Crystal.clear();
}
void TExlData::Dump() const
{
cout << "XXXXXXXXXXXXXXXXXXXXXXXX New Event XXXXXXXXXXXXXXXXX" << endl;
cout << "XXXXXXXXXXXXXXXXXXXXXXXX New Event XXXXXXXXXXXXXXXXX" << endl;
cout << "Energy:" << endl;
for(unsigned short i = 0 ; i<fExlE_Energy.size() ; i ++)
{
cout << "Exl Number " << fExlE_Number[i]
<< " Crystal Number " << fExlE_Crystal[i]
<< " Energy: " << fExlE_Energy[i] << endl;
}
cout << "Time:" << endl;
for(unsigned short i = 0 ; i<fExlT_Time.size() ; i ++)
{
cout << "Exl Number " << fExlT_Number[i]
<< " Crystal Number " << fExlT_Crystal[i]
<< " Time: " << fExlT_Time[i] << endl;
}
for(unsigned short i = 0 ; i<fExl_Energy.size() ; i ++)
{
cout << "Exl Number " << fExl_Number[i] << " Crystal Number " << fExl_Crystal[i] << " Energy: " << fExl_Energy[i] << " Time: " << fExl_Time[i] << endl;
}
}
double TExlData::GetEnergyCrystalNumber(const int crystal){
double energy = -1000;
unsigned int sizeE = GetEMult();
for (unsigned int i = 0; i < sizeE; i++){
if (GetE_CrystalNumber(i) == crystal){
energy = GetEnergy(i);
// return GetEnergy(i);
}
}
return energy;
}
double TExlData::GetTimeCrystalNumber(const int crystal){
double time = -1000;
unsigned int sizeT = GetTMult();
for (unsigned int i = 0; i < sizeT; i++){
if (GetT_CrystalNumber(i)==crystal){
time = GetTime(i);
// return GetTime(i);
}
}
return time;
}
}
......@@ -29,38 +29,58 @@ using namespace std ;
class TExlData : public TObject {
private:
// ADC
vector<double> fExl_Energy;
vector<double> fExl_Time;
vector<short> fExl_Number ;
vector<short> fExl_Crystal;
private:
// ADC
vector<double> fExlE_Energy;
vector<short> fExlE_Number ;
vector<short> fExlE_Crystal;
vector<double> fExlT_Time;
vector<short> fExlT_Number ;
vector<short> fExlT_Crystal;
public:
TExlData();
virtual ~TExlData();
void Clear();
void Clear(const Option_t*) {};
void Dump() const;
public:
TExlData();
virtual ~TExlData();
///////////////////// GETTERS ////////////////////////
inline double GetEnergy(const int& i) const { return fExl_Energy[i] ;}
inline double GetTime(const int& i) const { return fExl_Time[i] ;}
inline int GetExlNumber(const int& i) const { return fExl_Number[i] ;}
inline int GetCrystalNumber(const int& i) const { return fExl_Crystal[i] ;}
//Mult
inline double GetMult() const { return fExl_Energy.size() ;}
void Clear();
void Clear(const Option_t*) {};
void Dump() const;
///////////////////// SETTERS ////////////////////////
inline void SetEandTime(const int& N,const int& C,const double& E,const double& T){
fExl_Energy.push_back(E);
fExl_Time.push_back(T);
fExl_Number.push_back(N) ;
fExl_Crystal.push_back(C);
}
//////////////////////////////////////////////////////////////
// Getters and Setters
// Prefer inline declaration to avoid unnecessary called of
// frequently used methods
// add //! to avoid ROOT creating dictionnary for the methods
public:
///////////////////// GETTERS ////////////////////////
inline double GetEnergy(const int& i) const { return fExlE_Energy[i] ;}
double GetEnergyCrystalNumber(const int crystal);
double GetTimeCrystalNumber(const int crystal);
inline double GetTime(const int& i) const { return fExlT_Time[i] ;}
inline int GetE_ExlNumber(const int& i) const { return fExlE_Number[i] ;}
inline int GetE_CrystalNumber(const int& i) const { return fExlE_Crystal[i] ;}
inline int GetT_ExlNumber(const int& i) const { return fExlT_Number[i] ;}
inline int GetT_CrystalNumber(const int& i) const { return fExlT_Crystal[i] ;}
ClassDef(TExlData,1) // ExlData structure
//Mult
inline unsigned int GetEMult() const {return fExlE_Energy.size();}
inline unsigned int GetTMult() const {return fExlT_Time.size();}
///////////////////// SETTERS ////////////////////////
inline void SetEnergy(const int& N,const int& C,const double& E){
fExlE_Number.push_back(N) ;
fExlE_Crystal.push_back(C);
fExlE_Energy.push_back(E);
}
inline void SetTime(const int& N,const int& C,const double& T){
fExlT_Number.push_back(N) ;
fExlT_Crystal.push_back(C);
fExlT_Time.push_back(T);
}
ClassDef(TExlData,1) // ExlData structure
};
#endif
......@@ -7,21 +7,24 @@
/*****************************************************************************
* Original Author: L. Lefebvre contact address: lefebvrl@ipno.in2p3.fr *
* Update : V. Girard-Alcindor contact address: girardalcindor@ganil.fr *
* *
* Creation Date : October 2011 *
* Last update : *
* Last update : 2018 *
*---------------------------------------------------------------------------*
* Decription: *
* This class hold the Exl Detector Physics *
* *
*---------------------------------------------------------------------------*
* Comment: *
* Comment: Due to the way the detector is build, one has to take into *
* account the cross-talk between two neighboring crystals. *
* *
* *
*****************************************************************************/
// NPL
#include "TExlPhysics.h"
#include "TExlData.h"
#include "RootOutput.h"
#include "RootInput.h"
#include "NPDetectorFactory.h"
......@@ -38,199 +41,320 @@ using namespace std;
#include "TChain.h"
// tranform an integer to a string
string itoa(int value)
{
char buffer [33];
sprintf(buffer,"%d",value);
return buffer;
string itoa(int value){
char buffer [33];
sprintf(buffer,"%d",value);
return buffer;
}
ClassImp(TExlPhysics)
///////////////////////////////////////////////////////////////////////////
TExlPhysics::TExlPhysics()
{
NumberOfCrystal = 0;
EventData = new TExlData ;
EventPhysics = this ;
// Raw Threshold
m_E_RAW_Threshold = 0;
}
ClassImp(TExlPhysics);
///////////////////////////////////////////////////////////////////////////
TExlPhysics::~TExlPhysics()
{}
TExlPhysics::TExlPhysics(){
m_NumberOfDetector=0;
m_NumberOfCrystal = 18;
m_EventData = new TExlData;
m_EventPhysics = this;
m_E_RAW_Threshold = 0;
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::Clear()
{
CrystalNumber.clear() ;
EXL_Energy.clear() ;
EXL_Time=-1000 ;
}
A = 0.;
B = 0.;
a = 0.;
b = 0.;
v_A.clear();
v_B.clear();
v_a.clear();
v_b.clear();
m_Crystal_i = 0;
m_NeighborCrystal_i = 0;
m_Channel_i = 0.;
m_NeighborChannel_i = 0.;
m_Energy_i = 0.;
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::ReadConfiguration(NPL::InputParser parser) {
vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("EXL");
if(NPOptionManager::getInstance()->GetVerboseLevel())
cout << "//// " << blocks.size() << " detectors found " << endl;
vector<string> token= {"POS"};
vector<TVector3> Center_CsI_Crystals;
for(unsigned int i = 0 ; i < blocks.size() ; i++){
if(blocks[i]->HasTokenList(token)){
TVector3 pos = blocks[i]->GetTVector3("POS","mm");
Center_CsI_Crystals.push_back(pos);
}
// std::ifstream pre_calibration_coeff;
// pre_calibration_coeff.open("/home/alcindor/Documents/Projects/EXL/Analysis/Calibration/pre_calibration_coeff2.txt");
else{
cout << "ERROR: check your input file formatting " << endl;
exit(1);
// std::ifstream calibration_coeff;
// calibration_coeff.open("/home/alcindor/Documents/Projects/EXL/Analysis/Calibration/Co60_EXL_Calibration.txt");
// if(v_A.size()==0){
// while(pre_calibration_coeff >> B >> A){
// v_A.push_back(A);
// v_B.push_back(B);
// }
// }
// if(v_a.size()==0){
// while(calibration_coeff >> Channel_name >> b >> a){
// v_a.push_back(a);
// v_b.push_back(b);
// }
// }
for (int i = 0; i < 30 ; i++){
v_A.push_back(1);
v_B.push_back(1);
v_a.push_back(1);
v_b.push_back(1);
}
}
AddEXL(Center_CsI_Crystals) ;
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::AddEXL(vector <TVector3> Center_CsI_Crystals)
{
// In progress
// Needs the detectors dimensions
// Crystals in this detector can be mapped differently for every module
// Rows_y (0,1,2) Cols_x (0,1,2,3,4,5)
//local zero is at row 1, col 2.5
double arbit_length = 100; // detector arbit. face length and width, 100 mm
double col_pitch = arbit_length/6;
double row_pitch = arbit_length/3;
for(unsigned i=0; i<Center_CsI_Crystals.size(); i++)
{
TVector3 direction = Center_CsI_Crystals.at(i).Unit();
//place all 18 crystals in place
for(unsigned iCry=0; iCry<18; iCry++)
{
int Row_y = iCry/6;
int Col_x = iCry%6;
// Define Detector position localy
TVector3 localPos( (Row_y-1)*row_pitch,(Col_x-2.5)*row_pitch,0 );
// Rotate
localPos.RotateUz(direction);
TVector3 globalPos = Center_CsI_Crystals.at(i)+localPos;
cout << globalPos.X() << " " << globalPos.Y() << " " << globalPos.Z() << endl;
CsIPosition.push_back(globalPos);
}
}
TExlPhysics::~TExlPhysics(){
}
///////////////////////////////////////////////////////////////////////////
TVector3 TExlPhysics::GetPositionOfInteraction(int det, int cry)
{
int N = det*18+cry; // hyper crystal number
TVector3 Position = TVector3 (CsIPosition.at(N).X(),CsIPosition.at(N).Y(),CsIPosition.at(N).Z()) ;
void TExlPhysics::Clear(){
CrystalNumber.clear() ;
Energy.clear() ;
RawEnergy.clear();
Time.clear() ;
AddBack.clear();
}
return(Position) ;
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::ReadConfiguration(NPL::InputParser parser){
vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("EXL");
if(NPOptionManager::getInstance()->GetVerboseLevel())
cout << "//// " << blocks.size() << " detectors found " << endl;
}
vector<string> token = {"POS"};
vector<TVector3> Center_CsI_Crystals;
for(unsigned int i = 0; i < blocks.size(); i++){
if(blocks[i]->HasTokenList(token)){
TVector3 pos = blocks[i]->GetTVector3("POS","mm");
Center_CsI_Crystals.push_back(pos);
}
else{
cout << "ERROR: check your input file formatting " << endl;
exit(1);
}
}
AddEXL(Center_CsI_Crystals) ;
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::AddParameterToCalibrationManager()
{
CalibrationManager* Cal = CalibrationManager::getInstance();
void TExlPhysics::AddEXL(vector <TVector3> Center_CsI_Crystals){
m_NumberOfDetector++;
// In progress
// Needs the detectors dimensions
// Crystals in this detector can be mapped differently for every module
// Rows_y (0,1,2) Cols_x (0,1,2,3,4,5)
// local zero is at row 1, col 2.5
double length = 110; // detector arbit. face length and width, 100 mm
double width = 7.2;
double height = 9.8;
for(unsigned int i = 0 ; i < NumberOfCrystal ; i++)
{
Cal->AddParameter("EXL", "_E_"+ NPL::itoa(i+1),"EXL_E_"+ NPL::itoa(i+1)) ;
Cal->AddParameter("EXL", "_T_"+ NPL::itoa(i+1),"EXL_T_"+ NPL::itoa(i+1)) ;
}
}
for(unsigned i=0; i<Center_CsI_Crystals.size(); i++){
TVector3 direction = Center_CsI_Crystals.at(i).Unit();
//place all 18 crystals in place
for(unsigned iCry=0; iCry<18; iCry++){
int Row_y = iCry/6;
int Col_x = iCry%6;
// Define Detector position localy
TVector3 localPos((Row_y-1)*length,(Col_x-2.5)*length,0);
// Rotate
localPos.RotateUz(direction);
TVector3 globalPos = Center_CsI_Crystals.at(i)+localPos;
CsIPosition.push_back(globalPos);
}
}
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::InitializeRootInputRaw()
{
TChain* inputChain = RootInput::getInstance()->GetChain() ;
inputChain->SetBranchStatus ( "Exl" , true ) ;
inputChain->SetBranchStatus ( "fExl_*" , true ) ;
inputChain->SetBranchAddress( "Exl" , &EventData ) ;
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::InitializeRootInputPhysics()
{
TChain* inputChain = RootInput::getInstance()->GetChain();
inputChain->SetBranchStatus ( "Exl", true );
inputChain->SetBranchStatus ( "CrystalNumber", true );
inputChain->SetBranchStatus ( "EXL_Energy", true );
inputChain->SetBranchStatus ( "EXL_Time", true );
inputChain->SetBranchAddress( "Exl", &EventPhysics );
}
TVector3 TExlPhysics::GetPositionOfInteraction(int det, int cry){
int N = det*18+cry; // hyper crystal number
TVector3 Position = TVector3(CsIPosition.at(N).X(),
CsIPosition.at(N).Y(),
CsIPosition.at(N).Z());
return(Position) ;
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::InitializeRootOutput()
{
TTree* outputTree = RootOutput::getInstance()->GetTree() ;
outputTree->Branch( "Exl" , "TExlPhysics" , &EventPhysics ) ;
}
void TExlPhysics::AddParameterToCalibrationManager(){
CalibrationManager* Cal = CalibrationManager::getInstance();
for(unsigned int i = 0 ; i < m_NumberOfCrystal ; i++){
Cal->AddParameter("EXL", "E_"+ NPL::itoa(i+1), "EXL_E_" + NPL::itoa(i+1));
Cal->AddParameter("EXL", "T_"+ NPL::itoa(i+1), "EXL_T_" + NPL::itoa(i+1));
}
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::BuildPhysicalEvent()
{
BuildSimplePhysicalEvent() ;
}
void TExlPhysics::InitializeRootInputRaw(){
TChain* inputChain = RootInput::getInstance()->GetChain();
inputChain->SetBranchStatus("EXL", true);
inputChain->SetBranchStatus("fExl_*", true);
inputChain->SetBranchAddress("EXL", &m_EventData );
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::BuildSimplePhysicalEvent(){
if(EventData->GetMult()>0){
for(unsigned int i = 0 ; i < EventData->GetMult() ; i++){
if(EventData->GetEnergy(i)>m_E_RAW_Threshold){
CrystalNumber.push_back( EventData->GetCrystalNumber(i)) ;
EXL_Energy.push_back( CalibrationManager::getInstance()->ApplyCalibration("EXL/_E_" + NPL::itoa( EventData->GetExlNumber(i) ),EventData->GetEnergy(i) ) );
}
}
}
else{
CrystalNumber.push_back(-1000) ;
EXL_Energy.push_back(-1000);
}
if(EventData->GetMult()>0){
for(unsigned int i = 0 ; i < EventData->GetMult() ; i++){
//EXL_Time.push_back( CalibrationManager::getInstance()->ApplyCalibration("EXL/_T_" + NPL::itoa( EventData->GetExlNumber(i) ),EventData->GetTime(i) ) );
EXL_Time=EventData->GetTime(i);
}
}
else{
EXL_Time=-1000;
}
void TExlPhysics::InitializeRootInputPhysics(){
TChain* inputChain = RootInput::getInstance()->GetChain();
inputChain->SetBranchStatus("EXL", true);
inputChain->SetBranchStatus("CrystalNumber", true);
inputChain->SetBranchStatus("EXL_Energy", true);
inputChain->SetBranchStatus("EXL_AddBack", true );
inputChain->SetBranchStatus("EXL_Time", true );
inputChain->SetBranchAddress("EXL", &m_EventPhysics);
inputChain->SetBranchAddress("EXL", &m_EventData );
}
///////////////////////////////////////////////////////////////////////////
void TExlPhysics::InitializeRootOutput(){
TTree* outputTree = RootOutput::getInstance()->GetTree();
outputTree->Branch("EXL", "TExlPhysics", &m_EventPhysics);
}
///////////////////////////////////////////////////////////////////////////
double TExlPhysics::DopplerCorrection(double E, double Theta, double beta)
{
double Pi = 3.141592654 ;
void TExlPhysics::BuildPhysicalEvent(){
double E_corr = 0;
double gamma = 1./sqrt(1-beta*beta);
unsigned int sizeE = m_EventData->GetEMult();
for(unsigned int i = 0 ; i < sizeE ; i++){
E_corr = gamma*E*(1.-beta*cos(Theta*Pi/180.));
CrystalNumber.push_back(m_EventData->GetE_CrystalNumber(i));
return(E_corr);
double energy = CalibrationManager::getInstance()->ApplyCalibration(
"EXL/E_" + NPL::itoa(m_EventData->GetE_CrystalNumber(i)),
m_EventData->GetEnergy(i));
}
Energy.push_back(energy);
////////////////////////////////////////////////////////////////////////////////
// Construct Method to be pass to the DetectorFactory //
////////////////////////////////////////////////////////////////////////////////
NPL::VDetector* TExlPhysics::Construct(){
return (NPL::VDetector*) new TExlPhysics();
}
// look for associated time
bool check = false;
unsigned int sizeT = m_EventData->GetTMult();
for(unsigned int j = 0 ; j < sizeT ; j++){
if(m_EventData->GetT_CrystalNumber(j)==m_EventData->GetE_CrystalNumber(i)){
double time = CalibrationManager::getInstance()->ApplyCalibration(
"EXL/T_" + NPL::itoa(m_EventData->GetT_CrystalNumber(j)),
m_EventData->GetTime(j));
// cout << time << endl;
////////////////////////////////////////////////////////////////////////////////
// Registering the construct method to the factory //
////////////////////////////////////////////////////////////////////////////////
extern "C"{
class proxy_exl{
public:
proxy_exl(){
NPL::DetectorFactory::getInstance()->AddToken("EXL","Exl");
NPL::DetectorFactory::getInstance()->AddDetector("EXL",TExlPhysics::Construct);
Time.push_back(time);
check = true;
break;
}
}
};
if(!check)
Time.push_back(-1000);
proxy_exl p;
}
}
// for(int m_Crystal_i = 1; m_Crystal_i < 19; m_Crystal_i++){
// // Getting NeighborCrystalNumber
// if(m_Crystal_i % 2 == 0){
// m_NeighborCrystal_i = m_Crystal_i - 1;
// }
// else{
// m_NeighborCrystal_i = m_Crystal_i + 1;
// }
// // Getting channel of both crystals
// m_Channel_i = m_EventData->GetEnergyCrystalNumber(m_Crystal_i);
// m_NeighborChannel_i = m_EventData->GetEnergyCrystalNumber(m_NeighborCrystal_i);
// // Getting time
// Time.push_back(m_EventData->GetTimeCrystalNumber(m_Crystal_i));
// if (m_Channel_i != -10){
<