Skip to content
Snippets Groups Projects
Commit 13436990 authored by Hugo Jacob's avatar Hugo Jacob
Browse files

Modified CATS mask calibration code

parent 89330b3d
No related branches found
No related tags found
1 merge request!27Draft: [Epic] Preparation of the environement for the new GaseousDetectorScorers...
...@@ -223,9 +223,10 @@ void TCATSPhysics::BuildPhysicalEvent(){ ...@@ -223,9 +223,10 @@ void TCATSPhysics::BuildPhysicalEvent(){
// }); // });
// At least two CATS need to gave back position in order to reconstruct on Target // At least two CATS need to gave back position in order to reconstruct on Target
if(Positions.size()>1){ if(Positions.size()>1){
double t = (m_Zproj-Positions[2].first)/(Positions[2].first-Positions[1].first); double t = (m_Zproj-Positions[2].first)/(m_Zproj-Positions[1].first);
PositionOnTargetX= Positions[2].second.first + (Positions[2].second.first-Positions[1].second.first)*t; PositionOnTargetX= (Positions[2].second.first - Positions[1].second.first*t)/(1 - t);
PositionOnTargetY= Positions[2].second.second + (Positions[2].second.second-Positions[1].second.second)*t; PositionOnTargetY= (Positions[2].second.second - Positions[1].second.second*t)/(1 - t);
// PositionOnTargetY= Positions[2].second.second + (Positions[2].second.second-3-Positions[1].second.second)*t;
if(Mask1_Z != 0 && Mask2_Z != 0) if(Mask1_Z != 0 && Mask2_Z != 0)
{ {
double tmask1 = (Mask1_Z-Positions[2].first)/(Positions[2].first-Positions[1].first); double tmask1 = (Mask1_Z-Positions[2].first)/(Positions[2].first-Positions[1].first);
......
...@@ -75,13 +75,14 @@ void Analysis::Init(){ ...@@ -75,13 +75,14 @@ void Analysis::Init(){
bool Analysis::UnallocateBeforeBuild(){ bool Analysis::UnallocateBeforeBuild(){
// std::cout << "test unallocate" << std::endl; // std::cout << "test unallocate" << std::endl;
return true;
GATCONFMASTER = **GATCONFMASTER_; GATCONFMASTER = **GATCONFMASTER_;
return (GATCONFMASTER > 0); return (GATCONFMASTER > 0);
//return true; //return true;
} }
bool Analysis::UnallocateBeforeTreat(){ bool Analysis::UnallocateBeforeTreat(){
for(int i = 0; i < 10; i++){ /* for(int i = 0; i < 10; i++){
PlasticRaw[i] = (*PlasticRaw_ )[i]; PlasticRaw[i] = (*PlasticRaw_ )[i];
PlasticRawTS[i] = (*PlasticRaw_TS_)[i]; PlasticRawTS[i] = (*PlasticRaw_TS_)[i];
} }
...@@ -125,6 +126,7 @@ bool Analysis::UnallocateBeforeTreat(){ ...@@ -125,6 +126,7 @@ bool Analysis::UnallocateBeforeTreat(){
TAC_PL_4TS = **TAC_PL_4_TS_; TAC_PL_4TS = **TAC_PL_4_TS_;
TAC_PL_5 = **TAC_PL_5_; TAC_PL_5 = **TAC_PL_5_;
TAC_PL_5TS = **TAC_PL_5_TS_; TAC_PL_5TS = **TAC_PL_5_TS_;
*/
return true; return true;
} }
...@@ -137,8 +139,7 @@ bool Analysis::FillOutputCondition(){ ...@@ -137,8 +139,7 @@ bool Analysis::FillOutputCondition(){
void Analysis::TreatEvent(){ void Analysis::TreatEvent(){
// if(M2->CsI_E.size() > 0) // if(M2->CsI_E.size() > 0)
// std::cout << "Analysis test " << M2->CsI_E[0] << " " << M2->CsI_E.size() << " " << "\n \n"; // std::cout << "Analysis test " << M2->CsI_E[0] << " " << M2->CsI_E.size() << " " << "\n \n";*
ReInit(); ReInit();
// std::cout << CATS->PositionX.size() << std::endl; // std::cout << CATS->PositionX.size() << std::endl;
//////////////////// MUST2 Part //////////////////// //////////////////// MUST2 Part ////////////////////
...@@ -158,7 +159,7 @@ void Analysis::TreatEvent(){ ...@@ -158,7 +159,7 @@ void Analysis::TreatEvent(){
//BeamImpact = TVector3(0,0,0); //BeamImpact = TVector3(0,0,0);
BeamDirection = TVector3(CATS->PositionX[1] - CATS->PositionX[0],CATS->PositionY[1] - CATS->PositionY[0],CATS->PositionZ[1] - CATS->PositionZ[0]); BeamDirection = TVector3(CATS->PositionX[0] - CATS->PositionX[1],CATS->PositionY[0] - CATS->PositionY[1],CATS->PositionZ[0] - CATS->PositionZ[1]);
// std::cout << "Position XY " << CATS->PositionX[1] - CATS->PositionX[0] << " " << CATS->PositionY[1] - CATS->PositionY[0] << " " << CATS->PositionZ[1] - CATS->PositionZ[0] << std::endl; // std::cout << "Position XY " << CATS->PositionX[1] - CATS->PositionX[0] << " " << CATS->PositionY[1] - CATS->PositionY[0] << " " << CATS->PositionZ[1] - CATS->PositionZ[0] << std::endl;
// BeamDirection = TVector3(0,0,1); // BeamDirection = TVector3(0,0,1);
...@@ -200,7 +201,6 @@ void Analysis::TreatEvent(){ ...@@ -200,7 +201,6 @@ void Analysis::TreatEvent(){
Energy[ParticleType[i]] = LightTarget[ParticleType[i]].EvaluateInitialEnergy(Energy[ParticleType[i]] ,TargetThickness*0.5, ThetaNormalTarget); Energy[ParticleType[i]] = LightTarget[ParticleType[i]].EvaluateInitialEnergy(Energy[ParticleType[i]] ,TargetThickness*0.5, ThetaNormalTarget);
else else
Energy[ParticleType[i]] = -1000; Energy[ParticleType[i]] = -1000;
// std::cout << "DILO tata " << Energy[ParticleType[i]] << std::endl;
} }
...@@ -214,6 +214,7 @@ void Analysis::TreatEvent(){ ...@@ -214,6 +214,7 @@ void Analysis::TreatEvent(){
// Part 3 : Excitation Energy Calculation // Part 3 : Excitation Energy Calculation
M2_Ex_p.push_back(reaction->ReconstructRelativistic( Energy["proton"] , M2_ThetaLab[countMust2] )); M2_Ex_p.push_back(reaction->ReconstructRelativistic( Energy["proton"] , M2_ThetaLab[countMust2] ));
M2_Ex_d.push_back(Reaction_pd->ReconstructRelativistic( Energy["deuteron"] , M2_ThetaLab[countMust2] )); M2_Ex_d.push_back(Reaction_pd->ReconstructRelativistic( Energy["deuteron"] , M2_ThetaLab[countMust2] ));
std::cout << "oui " << M2_Ex_d[countMust2] << std::endl;
M2_Ex_t.push_back(Reaction_pt->ReconstructRelativistic( Energy["triton"] , M2_ThetaLab[countMust2] )); M2_Ex_t.push_back(Reaction_pt->ReconstructRelativistic( Energy["triton"] , M2_ThetaLab[countMust2] ));
M2_Ex_a.push_back(Reaction_p3He->ReconstructRelativistic( Energy["alpha"] , M2_ThetaLab[countMust2] )); M2_Ex_a.push_back(Reaction_p3He->ReconstructRelativistic( Energy["alpha"] , M2_ThetaLab[countMust2] ));
...@@ -232,6 +233,47 @@ void Analysis::TreatEvent(){ ...@@ -232,6 +233,47 @@ void Analysis::TreatEvent(){
// //
}//end loop MUST2 }//end loop MUST2
} }
/*
}
if(M2->Si_E.size() > 0){
if(Inner6MVM > 0){
int ExoMultMax = 2;
if(OutersVM < 2){
ExoMultMax == OutersVM;
}
for(int ExoMult = 0; ExoMult < ExoMultMax;ExoMult++){
ExogamDetNb[ExoMult] = (OutersVN[ExoMult] - 32)/16;
CristalNb[ExoMult] = (OutersVN[ExoMult] - (2+ ExogamDetNb[ExoMult])*16)/4;
SegmentNb[ExoMult] = (OutersVN[ExoMult] - 16*(ExogamDetNb[ExoMult] + 2) - 4*(CristalNb[ExoMult]));
}
if(Inner6MVM == 1){
if(OutersVM > 0 && BGOV[0] < 20){
Theta_seg = Exogam_Clovers_struc[ExogamDetNb[0]].Theta_Crystal_Seg[CristalNb[0]][SegmentNb[0]];
Phi_seg = Exogam_Clovers_struc[ExogamDetNb[0]].Phi_Crystal_Seg[CristalNb[0]][SegmentNb[0]];
EnergyDoppler = Doppler_Correction(Theta_seg,Phi_seg,0,0,Beta,Inner6MV[0]);
EnergyAddBackDoppler = EnergyDoppler;
}
}
if(Inner6MVM == 2 && OutersVM > 1 && BGOV[0] < 20 && BGOV[1] < 20){
if(Inner6MVN[0]/4 == Inner6MVN[1]/4){
EnergyAddBack = Inner6MV[0] + Inner6MV[1];
}
}
if(EnergyAddBack > -1000){
for(int i = 0;i < ExoMultMax; i++){
if(OutersV[i] > OutersV[highest_E]){
highest_E = i;
}
}
Theta_seg = Exogam_Clovers_struc[ExogamDetNb[highest_E]].Theta_Crystal_Seg[CristalNb[highest_E]][SegmentNb[highest_E]];
Phi_seg = Exogam_Clovers_struc[ExogamDetNb[highest_E]].Phi_Crystal_Seg[CristalNb[highest_E]][SegmentNb[highest_E]];
EnergyAddBackDoppler = Doppler_Correction(Theta_seg,Phi_seg,0,0,Beta,EnergyAddBack);
}
}
}
//for(unsigned int countMust2 = 0 ; countMust2 < M2->Si_E.size() ; countMust2++){ //for(unsigned int countMust2 = 0 ; countMust2 < M2->Si_E.size() ; countMust2++){
...@@ -254,11 +296,13 @@ void Analysis::TreatEvent(){ ...@@ -254,11 +296,13 @@ void Analysis::TreatEvent(){
// } // }
// } // }
//} //}
*/
} }
void Analysis::InitOutputBranch() { void Analysis::InitOutputBranch() {
/*
RootOutput::getInstance()->GetTree()->Branch("M2_TelescopeM",&M2_TelescopeM,"M2_TelescopeM/s"); RootOutput::getInstance()->GetTree()->Branch("M2_TelescopeM",&M2_TelescopeM,"M2_TelescopeM/s");
RootOutput::getInstance()->GetTree()->Branch("M2_CsI_E_p",&M2_CsI_E_p); RootOutput::getInstance()->GetTree()->Branch("M2_CsI_E_p",&M2_CsI_E_p);
RootOutput::getInstance()->GetTree()->Branch("M2_CsI_E_d",&M2_CsI_E_d); RootOutput::getInstance()->GetTree()->Branch("M2_CsI_E_d",&M2_CsI_E_d);
...@@ -278,7 +322,7 @@ void Analysis::InitOutputBranch() { ...@@ -278,7 +322,7 @@ void Analysis::InitOutputBranch() {
RootOutput::getInstance()->GetTree()->Branch("M2_Z",&M2_Z); RootOutput::getInstance()->GetTree()->Branch("M2_Z",&M2_Z);
RootOutput::getInstance()->GetTree()->Branch("M2_dE",&M2_dE); RootOutput::getInstance()->GetTree()->Branch("M2_dE",&M2_dE);
RootOutput::getInstance()->GetTree()->Branch("CsI_E_M2",&CsI_E_M2); RootOutput::getInstance()->GetTree()->Branch("CsI_E_M2",&CsI_E_M2);
RootOutput::getInstance()->GetTree()->Branch("M2_ECsI_from_deltaE",&M2_ECsI_from_deltaE); // RootOutput::getInstance()->GetTree()->Branch("M2_ECsI_from_deltaE",&M2_ECsI_from_deltaE);
RootOutput::getInstance()->GetTree()->Branch("GATCONF",&GATCONFMASTER); RootOutput::getInstance()->GetTree()->Branch("GATCONF",&GATCONFMASTER);
RootOutput:: getInstance()->GetTree()->Branch("TAC_CATS_PL",&TAC_CATS_PL,"TAC_CATS_PL/s"); RootOutput:: getInstance()->GetTree()->Branch("TAC_CATS_PL",&TAC_CATS_PL,"TAC_CATS_PL/s");
...@@ -321,15 +365,32 @@ void Analysis::InitOutputBranch() { ...@@ -321,15 +365,32 @@ void Analysis::InitOutputBranch() {
RootOutput:: getInstance()->GetTree()->Branch("IC_ZDDRaw",IC_ZDDRaw,"IC_ZDDRaw[6]/s"); RootOutput:: getInstance()->GetTree()->Branch("IC_ZDDRaw",IC_ZDDRaw,"IC_ZDDRaw[6]/s");
RootOutput:: getInstance()->GetTree()->Branch("IC_ZDDRawTS",IC_ZDDRawTS,"IC_ZDDRawTS[6]/l"); RootOutput:: getInstance()->GetTree()->Branch("IC_ZDDRawTS",IC_ZDDRawTS,"IC_ZDDRawTS[6]/l");
RootOutput:: getInstance()->GetTree()->Branch("EnergyDoppler",&EnergyDoppler,"EnergyDoppler/F");
RootOutput:: getInstance()->GetTree()->Branch("EnergyAddBack",&EnergyAddBack,"EnergyAddBack/F");
RootOutput:: getInstance()->GetTree()->Branch("EnergyAddBackDoppler",&EnergyAddBackDoppler,"EnergyAddBackDoppler/F");
RootOutput:: getInstance()->GetTree()->Branch("Inner6MVM",&Inner6MVM,"Inner6MVM/I");
RootOutput:: getInstance()->GetTree()->Branch("Inner6MV",Inner6MV,"Inner6MV[Inner6MVM]/F");
RootOutput:: getInstance()->GetTree()->Branch("Inner6MVN",Inner6MVN,"Inner6MVN[Inner6MVM]/s");
RootOutput:: getInstance()->GetTree()->Branch("Inner6MVTS",Inner6MVTS,"Inner6MVTS[Inner6MVM]/l");
RootOutput:: getInstance()->GetTree()->Branch("BGOVM",&BGOVM,"BGOVM/I");
RootOutput:: getInstance()->GetTree()->Branch("BGOV",BGOV,"BGOV[BGOVM]/F");
RootOutput:: getInstance()->GetTree()->Branch("BGOVN",BGOVN,"BGOVN[BGOVM]/s");
RootOutput:: getInstance()->GetTree()->Branch("DeltaTVM",&DeltaTVM,"DeltaTVM/I");
RootOutput:: getInstance()->GetTree()->Branch("DeltaTV",DeltaTV,"DeltaTV[DeltaTVM]/F");
RootOutput:: getInstance()->GetTree()->Branch("DeltaTVN",DeltaTVN,"DeltaTVN[DeltaTVM]/s");
RootOutput:: getInstance()->GetTree()->Branch("DeltaTVTS",DeltaTVTS,"DeltaTVTS[DeltaTVM]/l");
*/
} }
void Analysis::UnallocateVariables(){ void Analysis::UnallocateVariables(){
} }
void Analysis::InitInputBranch(){ void Analysis::InitInputBranch(){
TTreeReader* inputTreeReader = RootInput::getInstance()->GetTreeReader(); TTreeReader* inputTreeReader = RootInput::getInstance()->GetTreeReader();
GATCONFMASTER_ = new TTreeReaderValue<unsigned short>(*inputTreeReader,"GATCONFMASTER"); GATCONFMASTER_ = new TTreeReaderValue<unsigned short>(*inputTreeReader,"GATCONFMASTER");
//DATATRIG_CATS_ = new TTreeReaderValue<unsigned short>(*inputTreeReader,"DATATRIG_CATS"); /* //DATATRIG_CATS_ = new TTreeReaderValue<unsigned short>(*inputTreeReader,"DATATRIG_CATS");
PlasticRaw_ = new TTreeReaderArray<UShort_t>(*inputTreeReader,"PlasticRaw"); PlasticRaw_ = new TTreeReaderArray<UShort_t>(*inputTreeReader,"PlasticRaw");
PlasticRaw_TS_ = new TTreeReaderArray<ULong64_t>(*inputTreeReader,"PlasticRawTS"); PlasticRaw_TS_ = new TTreeReaderArray<ULong64_t>(*inputTreeReader,"PlasticRawTS");
...@@ -370,6 +431,20 @@ void Analysis::InitInputBranch(){ ...@@ -370,6 +431,20 @@ void Analysis::InitInputBranch(){
TAC_PL_4_TS_= new TTreeReaderValue<ULong64_t>(*inputTreeReader,"TAC_PL_4TS"); TAC_PL_4_TS_= new TTreeReaderValue<ULong64_t>(*inputTreeReader,"TAC_PL_4TS");
TAC_PL_5_= new TTreeReaderValue<UShort_t>(*inputTreeReader,"TAC_PL_5"); TAC_PL_5_= new TTreeReaderValue<UShort_t>(*inputTreeReader,"TAC_PL_5");
TAC_PL_5_TS_= new TTreeReaderValue<ULong64_t>(*inputTreeReader,"TAC_PL_5TS"); TAC_PL_5_TS_= new TTreeReaderValue<ULong64_t>(*inputTreeReader,"TAC_PL_5TS");
Inner6MVM_ = new TTreeReaderValue<int>(*inputTreeReader,"Inner6MRawM");
Inner6MV_ = new TTreeReaderArray<float>(*inputTreeReader,"Inner6MRaw");
Inner6MVN_ = new TTreeReaderArray<unsigned short>(*inputTreeReader,"Inner6MRawNr");
//RootInput:: getInstance()->GetChain()->SetBranchAddress("Inner6MVTS",Inner6MVTS);
OutersVM_ = new TTreeReaderValue<int>(*inputTreeReader,"OutersRawM");
OutersV_ = new TTreeReaderArray<float>(*inputTreeReader,"OutersRaw");
OutersVN_ = new TTreeReaderArray<unsigned short>(*inputTreeReader,"OutersRawNr");
BGOVM_ = new TTreeReaderValue<int>(*inputTreeReader,"BGORawM");
BGOV_ = new TTreeReaderArray<float>(*inputTreeReader,"BGORaw");
BGOVN_ = new TTreeReaderArray<unsigned short>(*inputTreeReader,"BGORawNr");
*/
} }
//////////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////////////
...@@ -392,7 +467,7 @@ void Analysis::ReInit(){ ...@@ -392,7 +467,7 @@ void Analysis::ReInit(){
M2_CsI_E_t.clear(); M2_CsI_E_t.clear();
M2_CsI_E_a.clear(); M2_CsI_E_a.clear();
M2_ECsI_from_deltaE.clear(); // M2_ECsI_from_deltaE.clear();
//ExNoBeam=ExNoProto.clear(); //ExNoBeam=ExNoProto.clear();
//EDC.clear(); //EDC.clear();
M2_ELab.clear(); M2_ELab.clear();
......
...@@ -105,9 +105,7 @@ class Analysis: public NPL::VAnalysis{ ...@@ -105,9 +105,7 @@ class Analysis: public NPL::VAnalysis{
std::vector<double> M2_Y; std::vector<double> M2_Y;
std::vector<double> M2_Z; std::vector<double> M2_Z;
std::vector<double> M2_dE; std::vector<double> M2_dE;
std::vector<double> M2_ECsI_from_deltaE ;
std::vector<double> Beta_from_deltaE;
std::vector<double> Beth_from_deltaE;
double OriginalBeamEnergy ; // AMEV double OriginalBeamEnergy ; // AMEV
double FinalBeamEnergy; double FinalBeamEnergy;
...@@ -133,7 +131,14 @@ class Analysis: public NPL::VAnalysis{ ...@@ -133,7 +131,14 @@ class Analysis: public NPL::VAnalysis{
unsigned short Inner6MVN[48]; unsigned short Inner6MVN[48];
TTreeReaderArray<unsigned short>* Inner6MVN_; TTreeReaderArray<unsigned short>* Inner6MVN_;
unsigned long long Inner6MVTS[48]; unsigned long long Inner6MVTS[48];
TTreeReaderArray<unsigned long long>* Inner6MVTS_;
int OutersVM;
TTreeReaderValue<int>* OutersVM_;
float OutersV[192];
TTreeReaderArray<float>* OutersV_;
unsigned short OutersVN[192];
TTreeReaderArray<unsigned short>* OutersVN_;
int BGOVM; int BGOVM;
TTreeReaderValue<int>* BGOVM_; TTreeReaderValue<int>* BGOVM_;
...@@ -169,19 +174,13 @@ class Analysis: public NPL::VAnalysis{ ...@@ -169,19 +174,13 @@ class Analysis: public NPL::VAnalysis{
float EnergyAddBackDoppler; float EnergyAddBackDoppler;
float EnergyAddBack; float EnergyAddBack;
int ExogamDetNb[3]; int ExogamDetNb[3];
// int CristalNb[3]; int CristalNb[3];
int SegmentNb[3]; int SegmentNb[3];
std::vector<int> event1; std::vector<int> event1;
std::vector<int> event2; std::vector<int> event2;
int highest_E; int highest_E;
int OutersVM;
TTreeReaderValue<int>* OutersVM_;
float OutersV[192];
TTreeReaderArray<float>* OutersV_;
unsigned short OutersVN[192];
TTreeReaderArray<unsigned short>* OutersVN_;
int DCRawM; int DCRawM;
unsigned short DCRaw[4]; unsigned short DCRaw[4];
......
...@@ -21,7 +21,7 @@ ...@@ -21,7 +21,7 @@
#include<iostream> #include<iostream>
using namespace std; using namespace std;
#include"Analysis.h" #include"Analysis_E805.h"
#include"NPAnalysisFactory.h" #include"NPAnalysisFactory.h"
#include"NPDetectorManager.h" #include"NPDetectorManager.h"
//////////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////////////
...@@ -218,8 +218,6 @@ bool Analysis::UnallocateBeforeBuild(){ ...@@ -218,8 +218,6 @@ bool Analysis::UnallocateBeforeBuild(){
//////////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////////////
void Analysis::TreatEvent(){ void Analysis::TreatEvent(){
}
/*
//std::cout << "Bonjour je suis Nicolas" << std::endl; //std::cout << "Bonjour je suis Nicolas" << std::endl;
ReInit(); ReInit();
...@@ -303,7 +301,7 @@ void Analysis::TreatEvent(){ ...@@ -303,7 +301,7 @@ void Analysis::TreatEvent(){
/////////////////////////////// Drift Chambers /////////////////// /////////////////////////////// Drift Chambers ///////////////////
/*if(DCRawM <=2 && DATATRIG_CATSTS[0] > 0){ if(DCRawM <=2 && DATATRIG_CATSTS[0] > 0){
for(int DCmult = 0; DCmult < DCRawM; DCmult++){ for(int DCmult = 0; DCmult < DCRawM; DCmult++){
DC_Nr = DCRawNr[DCmult]; DC_Nr = DCRawNr[DCmult];
...@@ -433,7 +431,6 @@ void Analysis::TreatEvent(){ ...@@ -433,7 +431,6 @@ void Analysis::TreatEvent(){
//std::cout << "Bonjour Nicolas est parti " << std::endl; //std::cout << "Bonjour Nicolas est parti " << std::endl;
} }
*/
......
CalibrationFilePath
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM1.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM2.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM3.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM4.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM1.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM2.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM3.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM4.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM1_r0183.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM1_r0184.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM2_r0183.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM2_r0184.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM3_r0189.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM3_r0190.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM4_r0187.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM4_r0188_manual.cal
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_alpha.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_alpha.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_alpha.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_alpha.txt
./calibration/CATS/CATS1X/CATSLin1X_r0164.txt
./calibration/CATS/CATS1Y/CATSLin1Y_r0164.txt
./calibration/CATS/CATS2X/CATSLin2X_r0166.txt
./calibration/CATS/CATS2Y/CATSLin2Y_r0166.txt
CalibrationFilePath
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM1.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM2.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM3.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_X_E_MM4.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM1.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM2.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM3.cal
./data/NPRoot/Calibration/CalAlpha_r0388_06+/Cal_Str_Y_E_MM4.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM1_r0183.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM1_r0184.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM2_r0183.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM2_r0184.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM3_r0189.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM3_r0190.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_X_T_MM4_r0187.cal
./data/NPRoot/Calibration/Time/Coeff/Cal_Str_Y_T_MM4_r0188_manual.cal
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_proton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_deuteron.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_triton.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM1_alpha.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM2_alpha.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM3_alpha.txt
./data/NPRoot/Calibration/Run_CSICalib/CsI_MM4_alpha.txt
./calibration/CATS/CATS1X/CATSLin1X_r0366.txt
./calibration/CATS/CATS1Y/CATSLin1Y_r0366.txt
./calibration/CATS/CATS2X/CATSLin2X_r0166.txt
./calibration/CATS/CATS2Y/CATSLin2Y_r0166.txt
...@@ -11,17 +11,17 @@ Target ...@@ -11,17 +11,17 @@ Target
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CATSDetector CATSDetector
CATSNumber = 1 CATSNumber = 1
X1_Y1= 43.16 -42.26 -1587.1 mm X1_Y1= 43.71 -41.76 -1587.1 mm
X28_Y1= -27.96 -42.26 -1587.1 mm X28_Y1= -27.41 -41.76 -1587.1 mm
X1_Y28= 43.16 28.86 -1587.1 mm X1_Y28= 43.71 29.36 -1587.1 mm
X28_Y28= -27.96 28.86 -1587.1 mm X28_Y28= -27.41 29.36 -1587.1 mm
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CATSDetector CATSDetector
CATSNumber = 2 CATSNumber = 2
X1_Y1= 33.66 -37.06 -1090.1 mm X1_Y1= 33.94 -37.06 -1090.1 mm
X28_Y1= -37.46 -37.06 -1090.1 mm X28_Y1= -37.18 -37.06 -1090.1 mm
X1_Y28= 33.66 34.06 -1090.1 mm X1_Y28= 33.94 34.06 -1090.1 mm
X28_Y28= -37.46 34.06 -1090.1 mm X28_Y28= -37.18 34.06 -1090.1 mm
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
MASK MASK
MaskNumber = 1 MaskNumber = 1
......
...@@ -7,7 +7,7 @@ Target ...@@ -7,7 +7,7 @@ Target
ANGLE= 0 deg ANGLE= 0 deg
X= 0 mm X= 0 mm
Y= 0 mm Y= 0 mm
Z= 10 mm Z= 0 mm
%%%%%%% Telescope 1 %%%%%%% %%%%%%% Telescope 1 %%%%%%%
M2Telescope M2Telescope
X1_Y1= -13.57 -104.78 299.83 mm X1_Y1= -13.57 -104.78 299.83 mm
...@@ -55,17 +55,17 @@ M2Telescope ...@@ -55,17 +55,17 @@ M2Telescope
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CATSDetector CATSDetector
CATSNumber = 1 CATSNumber = 1
X1_Y1= 43.16 -42.26 -1587.1 mm X1_Y1= 43.71 -41.76 -1587.1 mm
X28_Y1= -27.96 -42.26 -1587.1 mm X28_Y1= -27.41 -41.76 -1587.1 mm
X1_Y28= 43.16 28.86 -1587.1 mm X1_Y28= 43.71 29.36 -1587.1 mm
X28_Y28= -27.96 28.86 -1587.1 mm X28_Y28= -27.41 29.36 -1587.1 mm
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CATSDetector CATSDetector
CATSNumber = 2 CATSNumber = 2
X1_Y1= 33.66 -37.06 -1090.1 mm X1_Y1= 33.94 -37.06 -1090.1 mm
X28_Y1= -37.46 -37.06 -1090.1 mm X28_Y1= -37.18 -37.06 -1090.1 mm
X1_Y28= 33.66 34.06 -1090.1 mm X1_Y28= 33.94 34.06 -1090.1 mm
X28_Y28= -37.46 34.06 -1090.1 mm X28_Y28= -37.18 34.06 -1090.1 mm
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
MASK MASK
MaskNumber = 1 MaskNumber = 1
......
...@@ -17,8 +17,9 @@ using namespace std; ...@@ -17,8 +17,9 @@ using namespace std;
MCRMethod::MCRMethod() MCRMethod::MCRMethod()
{ {
Init(); //Init();
MetroHast(); // MetroHast();
Nappe();
//MinimizationFunction(); //MinimizationFunction();
} }
...@@ -185,8 +186,8 @@ void MCRMethod::Init(){ ...@@ -185,8 +186,8 @@ void MCRMethod::Init(){
Temperature = TemperatureInitiale; Temperature = TemperatureInitiale;
TemperatureFinale = 0.1; TemperatureFinale = 0.1;
Step = 0.1; Step = 0.1;
runmask[0] = "r0314_algo"; runmask[0] = "r0367_mask1";
runmask[1] = "r0315_algo"; runmask[1] = "r0368_mask1";
for(unsigned int i = 0; i < 2; i++){ for(unsigned int i = 0; i < 2; i++){
Chain[i] = new TChain("PhysicsTree"); Chain[i] = new TChain("PhysicsTree");
Chain[i]->Add(path+runmask[i]+".root"); //CATS1 Chain[i]->Add(path+runmask[i]+".root"); //CATS1
...@@ -281,8 +282,131 @@ double MCRMethod::MinimizationFunction(){ ...@@ -281,8 +282,131 @@ double MCRMethod::MinimizationFunction(){
return distdiff; return distdiff;
} }
void MCRMethod::Nappe(){
int X1 = 40;
int X2 = 40;
double stepx1 = 0.1;
double stepx2 = 0.1;
string Side = "Y";
TH1Map = new std::map<int,std::map<string,TH1F*>>;
double norm[2*X1+1][2*X2+1];
for(unsigned int i = 0; i < 2; i++){
for(int i1 = -X1; i1 <= X1; i1++){
for(int i2 = -X2; i2 <= X2; i2++){
(*TH1Map)[i][Form("h_m%i_x1%i_x2%i",i, i1, i2)] = new TH1F(Form("h_m%i_x1%i_x2%i",i, i1, i2),Form("h_m%i_x1%i_x2%i",i, i1, i2),200,-10,10);
}
}
}
// std::cout << "test2" << std::endl;
runmask[0] = "r0314_mask1";
runmask[1] = "r0315_mask1";
for(unsigned int i = 0; i < 2; i++){
Chain[i] = new TChain("PhysicsTree");
Chain[i]->Add(path+runmask[i]+".root"); //CATS1
if (!(Chain[i]->GetEntries())){
printf("ERROR in the Chain !! \n");
return;
}
// std::cout << "test3" << std::endl;
TreeReader[i] = new TTreeReader(Chain[i]);
CATSPhysics_[i]= new TTreeReaderValue<TCATSPhysics>(*TreeReader[i],"CATS");
}
TString CName = "CUT_314_mask1";
CFile[0] = new TFile(CName+".root");
CUT[0] = (TCutG*)CFile[0]->FindObjectAny(CName);
std::cout << CUT[0] << std::endl;
CName = "CUT_315_mask1";
CFile[1] = new TFile(CName+".root");
CUT[1] = (TCutG*)CFile[1]->FindObjectAny(CName);
std::cout << CUT[1] << std::endl;
for(unsigned int i = 0; i < 2; i++){
while (TreeReader[i]->Next()){
CATSPhysics = &**(CATSPhysics_[i]);
if(CATSPhysics->PositionX.size() == 2){
// std::cout << "test5 " << CATSPhysics->GetCATSMult() << " " << CATSPhysics->DetNumber[0] << " " << (CATSPhysics->PositionX).size() << " " << CATSPhysics->PositionY[i] << " " << CUT[i]->IsInside(CATSPhysics->PositionX[i],CATSPhysics->PositionY[i]) << std::endl;
if(CATSPhysics->DetNumber[0] == 1 && CUT[i]->IsInside(CATSPhysics->PositionX[i],CATSPhysics->PositionY[i])){
for(int i1 = -X1; i1 <= X1; i1++){
for(int i2 = -X2; i2 <= X2; i2++){
// std::cout << "test " << i1 << " " << i2 << std::endl;
// std::cout << "test4" << std::endl;
double x1, x2;
if(Side == "X"){
x1 = CATSPhysics->PositionX[0] +i1*stepx1;
x2 = CATSPhysics->PositionX[1] +i2*stepx2;
}
else if(Side == "Y"){
x1 = CATSPhysics->PositionY[0] +i1*stepx1;
x2 = CATSPhysics->PositionY[1] +i2*stepx2;
}
(*TH1Map)[i][Form("h_m%i_x1%i_x2%i",i, i1, i2)]->Fill(ProjectOnCats(i,x1,x2));
norm[X1+i1][X2+i2] = 0;
}
}
}
}
}
}
(*TH1Map)[0][Form("h_m%i_x1%i_x2%i",0, 0, 0)]->Draw();
(*TH1Map)[0][Form("h_m%i_x1%i_x2%i",0, 1, 1)]->Draw("same");
(*TH1Map)[0][Form("h_m%i_x1%i_x2%i",0, 2, 2)]->Draw("same");
auto c1 = new TCanvas;
TH2F* HeatMap[2];
auto fitfunc = new TF1("fitfunc","gausn",-10,10);
auto SumHeatMap = new TH2F("Sumheatmap","Sumheatmap",2*X1+1, -X1*stepx1, (X1+1)*stepx1, 2*X2+1, -X2*stepx2, (X2+1)*stepx2);
HeatMap[0] = new TH2F("heatmap_0","heatmap_0",2*X1+1, -X1*stepx1, (X1+1)*stepx1, 2*X2+1, -X2*stepx2, (X2+1)*stepx2);
HeatMap[1] = new TH2F("heatmap_1","heatmap_1",2*X1+1, -X1*stepx1, (X1+1)*stepx1, 2*X2+1, -X2*stepx2, (X2+1)*stepx2);
TSpectrum* Spec = new TSpectrum(1,1.);
for(unsigned int i = 0; i < 2; i++){
for(int i1 = -X1; i1 <= X1; i1++){
for(int i2 = -X2; i2 <= X2; i2++){
Int_t nfound = Spec->Search((*TH1Map)[i][Form("h_m%i_x1%i_x2%i",i, i1, i2)],2,"",0.15);
Double_t* PeaksPosX = Spec->GetPositionX();
fitfunc->SetParameters(200,PeaksPosX[0],1);
fitfunc->SetParLimits(1,-10,10);
(*TH1Map)[i][Form("h_m%i_x1%i_x2%i",i, i1, i2)]->Fit(fitfunc,"QL","",-10,10);
// std::cout << i1 << " " << i2 << " " << PeaksPosX[0] << std::endl;
double value;
if(Side == "X")
value = abs(fitfunc->GetParameter(1) - MaskPosX[i]);
else if(Side == "Y")
value = abs(fitfunc->GetParameter(1) - MaskPosY[i]);
HeatMap[i]->SetBinContent(X1+i1+1,X2+i2+1,value);
//HeatMap[i]->SetBinContent(i1+1,i2+1,value);
std::cout << norm[X1+i1][X2+i2] << std::endl;
norm[X1+i1][X2+i2] += value*value;
if(i ==1){
norm[X1+i1][X2+i2] = sqrt(norm[X1+i1][X2+i2]);
}
}
}
}
for(int i1 = -X1; i1 <= X1; i1++){
for(int i2 = -X2; i2 <= X2; i2++){
SumHeatMap->SetBinContent(X1+i1+1,X2+i2+1,norm[X1+i1][X2+i2]);
}
}
// SumHeatMap->Add(HeatMap[0],HeatMap[1],1,1);
for(unsigned int i = 0; i < 2; i++){
auto c = new TCanvas();
HeatMap[i]->Draw("colz");
}
auto c = new TCanvas();
SumHeatMap->Draw("colz");
}
float MCRMethod::ProjectOnCats(unsigned int i, double x1, double x2){
double tmask = (CATSPosZ[1]-MASKPosZ[i])/(CATSPosZ[0] - MASKPosZ[i]);
return x2 -x1*tmask;
};
......
...@@ -12,6 +12,9 @@ ...@@ -12,6 +12,9 @@
#include "TTreeReaderValue.h" #include "TTreeReaderValue.h"
#include "TH2F.h" #include "TH2F.h"
#include "TSpectrum2.h" #include "TSpectrum2.h"
#include "TSpectrum.h"
#include "TCutG.h"
#include "TF1.h"
// NPTOOL // NPTOOL
#include "TCATSPhysics.h" #include "TCATSPhysics.h"
...@@ -32,6 +35,10 @@ class MCRMethod ...@@ -32,6 +35,10 @@ class MCRMethod
void RandomStep(); void RandomStep();
double MinimizationFunction(); double MinimizationFunction();
void Init(); void Init();
float ProjectOnCats(unsigned int i, double x1, double x2);
void Nappe();
private: private:
...@@ -55,9 +62,11 @@ class MCRMethod ...@@ -55,9 +62,11 @@ class MCRMethod
private: private:
TChain *Chain[2]; TChain *Chain[2];
TTreeReader *TreeReader[2]; TTreeReader *TreeReader[2];
TString path = "./NPRootA/"; TString path = "./RootR/";
TTreeReaderValue<TCATSPhysics> *CATSPhysics_[2]; TTreeReaderValue<TCATSPhysics> *CATSPhysics_[2];
TCATSPhysics *CATSPhysics; TCATSPhysics *CATSPhysics;
TCutG *CUT[2];
TFile *CFile[2];
Double_t TopLeft[2][2]; Double_t TopLeft[2][2];
Double_t BotLeft[2][2]; Double_t BotLeft[2][2];
Double_t BotRight[2][2]; Double_t BotRight[2][2];
...@@ -66,8 +75,12 @@ class MCRMethod ...@@ -66,8 +75,12 @@ class MCRMethod
Double_t BotLeftGeo[4] = {-12.5,-12.3,-12.4,-13.6}; Double_t BotLeftGeo[4] = {-12.5,-12.3,-12.4,-13.6};
Double_t BotRightGeo[4] = {12.5,-12.3,12.6,-13.6}; Double_t BotRightGeo[4] = {12.5,-12.3,12.6,-13.6};
std::map<int,std::map<string,TH1F*>>* TH1Map;
TGraph* Graph; TGraph* Graph;
double MaskPosX[2] = {0,0.1};
double MaskPosY[2] = {0.2,-1.1};
double tini; double tini;
double tfinale; double tfinale;
double pas; double pas;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment