diff --git a/NPAnalysis/Gaspard/CrossSection.c b/NPAnalysis/Gaspard/CrossSection.c deleted file mode 100644 index 5640efb446f5464a98882b5a40ea015ef2745171..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/CrossSection.c +++ /dev/null @@ -1,121 +0,0 @@ -{ - gROOT->Reset(); - gStyle->SetOptStat(1); - gStyle->SetPalette(51,0); - /////////////////////// - /////////////////////// - - double Pi = 3.141592654 ; - -///Usefull myAngleInRad * RadToDeg = myAngle In deg :p - double DegToRad = Pi/180. ; // 2Pi/360 = Pi/180 - double RadToDeg = 180./Pi ; // 360/2Pi = 180/Pi - -TFile *file0 = TFile::Open("./Result/myResult.root"); - - cEA = new TCanvas("cEA","Kinematic Line" ,100,100,900,900); - hEA->Draw("COLZ"); - cEx = new TCanvas("cEx","Excitation Energy" ,100,100,600,600); - hEx->Draw(); - - cEHexa = new TCanvas("cEHexa","Hexaneutron bound Energy" ,100,100,600,600); - hEHexa->Draw(); - - cCM = new TCanvas("cCm" , "Cross Section (CM)" , 100 , 100 , 900, 900) ; - hThetaCM->Draw(); - - vector<double> relative_error; - double NumberOfHitX=0; - ifstream efficFile; - efficFile.open("./17cm.efficiency"); - string EffString, ThetaString; - double xxx = hThetaLabCM->GetNbinsY(); - - for(Int_t i = 1 ; i < hThetaLabCM->GetNbinsX() ; i++) - { - // Calculate Relative error - NumberOfHitX=hThetaLabCM->Integral(i, i, 1, xxx ); - if(NumberOfHitX!=0)relative_error.push_back( 1. / ( sqrt(NumberOfHitX) ) ); - else {relative_error.push_back(0) ;} - } - - for(Int_t j = 1 ; j < hThetaLabCM->GetNbinsY() ; j++) - { - //Efficacity correction - efficFile >> ThetaString >> EffString; - double Eff = atoi(EffString.c_str()); - for(Int_t i = 1 ; i < hThetaLabCM->GetNbinsX() ; i++) - { - if (Eff!=0) - { - double pipo = (hThetaLabCM->GetCellContent(i,j))/(Eff) ; - hThetaLabCM->SetCellContent( i, j, pipo ); - } - else hThetaLabCM->SetCellContent( i, j, 0 ); - } - } - - TH1D* hCrossSection = new TH1D(); - hCrossSection = hThetaLabCM->ProjectionX("Cross Section", 0, -1, "") ; - - for(Int_t i = 1 ; i < hCrossSection->GetNbinsX() ; i++) - { - // Calculate Error - hCrossSection->SetBinContent(i, - NumberOfHitX=hCrossSection->GetBinContent(i)/(sin(hCrossSection->GetBinCenter(i)*DegToRad) )); - double error=NumberOfHitX*relative_error[i]; - hCrossSection->SetBinError(i,error); - } - - hCrossSection->Draw(); - - string Path="../Data/CrossSection/11Li(d,3He)10He.txt"; - ifstream CrossSectionFile; - CrossSectionFile.open(Path.c_str()); - if( CrossSectionFile.is_open() ) - cout << " file " << Path << " loading " << endl; - else{ - cout << " Error, no file " << Path << " found" << endl; return;} - - - //Normalisation: - //Int_t Maximum_Bin = hCrossSection->GetMaximumBin() ; - Int_t Maximum_Bin = 3 ; - Double_t Maximum_theta = hCrossSection->GetBinCenter(Maximum_Bin) ; - Double_t Bin_Width = hCrossSection->GetBinWidth(Maximum_Bin) ; - Double_t Maximum = hCrossSection->GetBinContent(Maximum_Bin) ; - - Double_t Normalisation_Factor = 0 ; - Double_t ThetaVal = 0 ; - Double_t CSVal = 0 ; - TMarker marker ; - string theta,CS ; - - while( !CrossSectionFile.eof() && ThetaVal<60) - { - CrossSectionFile >> theta >> CS ; - - ThetaVal=atof(theta.c_str()) ; - - if(ThetaVal>(Maximum_theta-Bin_Width/1000) && ThetaVal<(Maximum_theta+Bin_Width/1000) ) - { - CSVal=atof(CS.c_str()) ; - Normalisation_Factor = Maximum/CSVal; - } - } - CrossSectionFile.close(); - CrossSectionFile.open(Path.c_str()); - - //Reading Cross Section - - ThetaVal=0; - while( !CrossSectionFile.eof() && ThetaVal<60) - { - - CrossSectionFile >> theta >> CS ; - CSVal=atof(CS.c_str()) ; - ThetaVal=atof(theta.c_str()) ; - marker->DrawMarker(ThetaVal,CSVal*Normalisation_Factor ) ; - } - -} diff --git a/NPAnalysis/Gaspard/Makefile b/NPAnalysis/Gaspard/Makefile new file mode 100644 index 0000000000000000000000000000000000000000..f5bfeb69ffd25b28c9ccb8a81eebb1b66664dfa2 --- /dev/null +++ b/NPAnalysis/Gaspard/Makefile @@ -0,0 +1,6 @@ + +Analyse: + make -C ./src + +clean: + make clean -C ./src diff --git a/NPAnalysis/Gaspard/Result/myResult.root b/NPAnalysis/Gaspard/Result/myResult.root deleted file mode 100644 index eb7aa8c9045793132a6c1c12c390e3b06e31eeb7..0000000000000000000000000000000000000000 Binary files a/NPAnalysis/Gaspard/Result/myResult.root and /dev/null differ diff --git a/NPAnalysis/Gaspard/Result/myResult_newGPD.root b/NPAnalysis/Gaspard/Result/myResult_newGPD.root deleted file mode 100644 index ab4d9b0461ed828ac52fcd4fab9f366248c0950e..0000000000000000000000000000000000000000 Binary files a/NPAnalysis/Gaspard/Result/myResult_newGPD.root and /dev/null differ diff --git a/NPAnalysis/Gaspard/RunToTreat.txt b/NPAnalysis/Gaspard/RunToTreat.txt new file mode 100644 index 0000000000000000000000000000000000000000..008354fd1b9fa2d1c15477a463be1e2bf76a6b7e --- /dev/null +++ b/NPAnalysis/Gaspard/RunToTreat.txt @@ -0,0 +1,4 @@ +TTreeName + SimulatedTree +RootFileName + ../../Outputs/Simulation/mySimul.root diff --git a/NPAnalysis/Gaspard/TimeOfFlight.c b/NPAnalysis/Gaspard/TimeOfFlight.c deleted file mode 100644 index 916932fc04c540832f14c0c4cf14b4889d8972be..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/TimeOfFlight.c +++ /dev/null @@ -1,14 +0,0 @@ -{ - - gROOT->Reset(); - gStyle->SetOptStat(1); - gStyle->SetPalette(51,0); - //gStyle->SetPalette(1); - /////////////////////// - /////////////////////// -TFile *file0 = TFile::Open("./Result/myResult.root"); - - TCanvas* cTOF = new TCanvas("cTOF","Time of Flight" ,10,10,800,800); - hTOF1234->Draw("COLZ"); - -} diff --git a/NPAnalysis/Gaspard/affich.c b/NPAnalysis/Gaspard/affich.c deleted file mode 100644 index 2c8537ea62be3b91ec797d490b630dd8f4ecc6eb..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/affich.c +++ /dev/null @@ -1,291 +0,0 @@ -{ - - gROOT->Reset(); - gStyle->SetOptStat(1); - //gStyle->SetPalette(51,0); - gStyle->SetPalette(1); - /////////////////////// - /////////////////////// -TFile *file0 = TFile::Open("./Result/myResult.root"); - -/* TCanvas* cTOF = new TCanvas("cTOF","Time of Flight" ,10,10,800,800); - cTOF->Divide(3,2); - - cTOF->cd(1); - hTOF1->Draw("COLZ"); - - cTOF->cd(2); - hTOF2->Draw("COLZ"); - - cTOF->cd(3); - hTOF3->Draw("COLZ"); - - cTOF->cd(4); - hTOF4->Draw("COLZ"); - - cTOF->cd(5); - hTOF5->Draw("COLZ"); */ - - double Pi = 3.141592654 ; - double DegToRad = Pi/180. ; // 2Pi/360 = Pi/180 - -/* vector<double> relative_error; - double NumberOfHitX=0; - cTheta2D = new TCanvas("cTheta2d","Cross Section" ,100,100,900,900); - hThetaLabCM->Draw("COLZ");*/ - cTheta = new TCanvas("cTheta","Cross Section" ,100,100,900,900); - hTheta->Draw(); -/* ifstream efficFile; - - efficFile.open("./efficiency/20cm2.efficiency"); - string EffString, ThetaString; - double xxx = hThetaLabCM->GetNbinsY(); - - for(Int_t i = 1 ; i < hThetaLabCM->GetNbinsX() ; i++) - { - // Calculate Relative error - NumberOfHitX=hThetaLabCM->Integral(i, i, 1, xxx ); - if(NumberOfHitX!=0)relative_error.push_back( ( sqrt(NumberOfHitX) ) / (NumberOfHitX) ); - else {relative_error.push_back(0) ;} - } - - for(Int_t j = 1 ; j < hThetaLabCM->GetNbinsY() ; j++) - { - //Efficacity correction - efficFile >> ThetaString >> EffString; - double Eff = atoi(EffString.c_str()); - for(Int_t i = 1 ; i < hThetaLabCM->GetNbinsX() ; i++) - { - if (Eff!=0) - { - double pipo = (hThetaLabCM->GetCellContent(i,j))/(Eff/100) ; - hThetaLabCM->SetCellContent( i, j, pipo ); - } - else hThetaLabCM->SetCellContent( i, j, 0 ); - } - } - - - TH1D* hCrossSection = new TH1D(); - hCrossSection = hThetaLabCM->ProjectionX("Cross Section", 0, -1, "") ; - - - for(Int_t i = 1 ; i < hCrossSection->GetNbinsX() ; i++) - { - // Calculate Error - NumberOfHitX=hCrossSection->GetBinContent(i); - double error=NumberOfHitX*relative_error[i]; - hCrossSection->SetBinError(i,error); - } - - hCrossSection->Draw(); - -/* ofstream efficFile; - efficFile.open("20cm.efficiency"); - - for(Int_t i = 1 ; i < hTheta->GetNbinsX() ; i++) - { - - Double_t Flux = 100000*sin(hTheta->GetBinCenter(i)*DegToRad)*2*DegToRad/2; - - if(Flux!=0) - { - Double_t Efficiency = ( ( Flux - hTheta->GetBinContent(i) ) / Flux ) * 100 ; - - hTheta->SetBinContent(i, 100-Efficiency ) ; - - efficFile << hTheta->GetBinCenter(i) << " " << 100-Efficiency << endl; - } - - else hTheta->SetBinContent(i, 0 ); - } - - - //hTheta->Draw();*/ - - //efficiency correction and error bar calculation -/* ifstream efficFile; - efficFile.open("./efficiency/20cmCM.efficiency"); - string EffString, ThetaString; - double Eff, relative_error, error, toto ; - - for(Int_t i = 1 ; i < hThetaCM->GetNbinsX() ; i++) - { - - efficFile >> ThetaString >> EffString; - Eff = atoi(EffString.c_str()); - if(hThetaCM->GetBinContent(i)!=0)relative_error = sqrt(hThetaCM->GetBinContent(i))/(hThetaCM->GetBinContent(i)); - - toto=hThetaCM->GetBinContent(i)*(100-Eff) ; - hThetaCM->SetBinContent(i, toto) ; - - error=hThetaCM->GetBinContent(i)*relative_error; - hThetaCM->SetBinError(i,error); - }*/ - - - - - //hThetaCM->Draw(); - -/* string Path="/home/Adrien/Desktop/geant/Simulation/Data/CrossSection/11Li(d,3He)10He.txt"; - ifstream CrossSectionFile; - CrossSectionFile.open(Path.c_str()); - if( CrossSectionFile.is_open() ) - cout << " file " << Path << " loading " << endl; - else{ - cout << " Error, no file " << Path << " found" << endl; return;} - - - //Normalisation: - Int_t Maximum_Bin = hCrossSection->GetMaximumBin() ; - Double_t Maximum_theta = hCrossSection->GetBinCenter(Maximum_Bin) ; - Double_t Bin_Width = hCrossSection->GetBinWidth(Maximum_Bin) ; - Double_t Maximum = hCrossSection->GetBinContent(Maximum_Bin) ; - - Double_t Normalisation_Factor=0; - - - - Double_t ThetaVal=0 ; - Double_t CSVal=0 ; - TMarker marker ; - string theta,CS ; - - while( !CrossSectionFile.eof() && ThetaVal<60) - { - CrossSectionFile >> theta >> CS ; - - ThetaVal=atof(theta.c_str()) ; - - if(ThetaVal>(Maximum_theta-Bin_Width/4) && ThetaVal<(Maximum_theta+Bin_Width/4) ) - { - CSVal=atof(CS.c_str()) ; - Normalisation_Factor = Maximum/CSVal; - } - } - CrossSectionFile.close(); - CrossSectionFile.open(Path.c_str()); - - //Reading Cross Section - - ThetaVal=0; - while( !CrossSectionFile.eof() && ThetaVal<60) - { - - CrossSectionFile >> theta >> CS ; - CSVal=atof(CS.c_str()) ; - ThetaVal=atof(theta.c_str()) ; - marker->DrawMarker(ThetaVal,CSVal*Normalisation_Factor) ; - } - - - -/* - //Normalisation: - Double_t Maximum_Bin = hTheta->GetMaximumBin() ; - Double_t Maximum = hTheta->GetBinContent(Maximum_Bin) ; - - Double_t Normalisation_Factor = Maximum/3.; - - - - //Reading Cross Section - Double_t ThetaVal=0 ; - Double_t CSVal=0 ; - TMarker marker ; - string theta,CS ; - - while( !CrossSectionFile.eof() && ThetaVal<60) - { - - CrossSectionFile >> theta >> CS ; - CSVal=atof(CS.c_str()) ; - ThetaVal=atof(theta.c_str()) ; - - ThetaVal=0.000479597*pow(ThetaVal,3)-0.0670771*pow(ThetaVal,2)+3.08149*ThetaVal-0.339958 ; - - marker->DrawMarker(ThetaVal,CSVal*Normalisation_Factor) ; - } -// cTheta->SetLogy(1); - - - - - /* - cKine = new TCanvas("cKine","Kinematics lines" ,10,10,800,800); - hKine->Draw("COLZ");*/ - -/* cKineInit = new TCanvas("cKineInit","Initial Kinematics lines" ,100,100,600,600); - hKineInit->Draw("COLZ"); - - cEDE = new TCanvas("EDE","EDE add Strip" ,100,100,600,600); - hEDEAddStrip->Draw("COLZ"); -/* cG = new TCanvas("cG","Strip position",500,100,800,600); - cG->Divide(3,2); - - cG->cd(1); - Agraph2D->Draw("P0") ; - cG->cd(2); - Bgraph2D->Draw("P0") ; - cG->cd(3); - Cgraph2D->Draw("P0") ; - cG->cd(4); - Dgraph2D->Draw("P0") ; - cG->cd(5); - Egraph2D->Draw("P0") ; - - cG2 = new TCanvas("cG2","all Strip position",500,100,800,600); - TOTgraph2D->Draw("P0") ; - */ - cH = new TCanvas("cH","Hit density",500,100,1000,800); - cH->Divide(3,2); - - cH->cd(1); - hHIT4->Draw("COLZ"); - cH->cd(2); - hHIT2->Draw("COLZ"); - cH->cd(3); - hHIT5->Draw("COLZ"); - cH->cd(4); - hHIT1->Draw("COLZ"); - cH->cd(5); - hHIT3->Draw("COLZ"); - cH->cd(6); - hThetaHeavy->Draw(); - -/* cH->cd(5); - hXZ->Draw("COLZ"); - cH->cd(6); - hXY->Draw("COLZ"); - - cEDE = new TCanvas("cEDE","EDE indentification Spectra",500,100,1000,800); - cEDE->Divide(3,2); - - cEDE->cd(1); - hEDE1->Draw("COLZ"); - cEDE->cd(2); - hEDE2->Draw("COLZ"); - cEDE->cd(3); - hEDE3->Draw("COLZ"); - cEDE->cd(4); - hEDE4->Draw("COLZ"); - cEDE->cd(5); - hEDE5->Draw("COLZ"); - cEDE->cd(6); - TH2F* hEDET= new TH2F("hEDET","",4000,-1,600, 400, -1, 25) ; - hEDET->Add(hEDE1); - hEDET->Add(hEDE2); - hEDET->Add(hEDE3); - hEDET->Add(hEDE4); - hEDET->Draw("CLOZ");*/ - - cEx = new TCanvas("cEx","Excitation Energy" ,100,100,600,600); - hEx->Draw(); - /*cE = new TCanvas("cE","Light Energy" ,500,100,1000,800); - cE->Divide(2); - cE->cd(1); - hE1234->Draw(); - cE->cd(2); - hE5->Draw();*/ -} diff --git a/NPAnalysis/Gaspard/efficiency.c b/NPAnalysis/Gaspard/efficiency.c deleted file mode 100644 index 9a04af185f75853d0bce5bae2009cd8cd01cae1a..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/efficiency.c +++ /dev/null @@ -1,58 +0,0 @@ -{ - gROOT->Reset(); - gStyle->SetOptStat(1); - gStyle->SetPalette(51,0); - /////////////////////// - /////////////////////// -TFile *file0 = TFile::Open("./Result/myResult.root"); - - double Pi = 3.141592654 ; - double DegToRad = Pi/180. ; // 2Pi/360 = Pi/180 - - - ofstream efficFile; - efficFile.open("12cm.efficiency"); - - for(Int_t i = 1 ; i < hTheta->GetNbinsX() ; i++) - { - Double_t Flux = 1000000./2. * sin(hTheta->GetBinCenter(i)*DegToRad) * (2*DegToRad) ; - - if(Flux!=0) - { - Double_t Efficiency = 100*( (hTheta->GetBinContent(i)) / Flux ); - hTheta->SetBinContent(i, Efficiency ) ; - } - - else hTheta->SetBinContent(i, 0 ); - - efficFile << hTheta->GetBinCenter(i) << " " << Efficiency << endl; - } - - cEff = new TCanvas("cEff","Efficiency" ,100,100,600,600); - hTheta->Draw(); - - cHit = new TCanvas("cHit","Hit" ,100,100,600,600); - hXY->Draw("COLZ"); - - cHit2 = new TCanvas("cHit2","Hit" ,100,100,600,600); - cHit2->Divide(3,2); - - cHit2->cd(1); - hXY1->Draw("COLZ"); - - cHit2->cd(2); - hXY2->Draw("COLZ"); - - cHit2->cd(3); - hXY3->Draw("COLZ"); - - cHit2->cd(4); - hXY4->Draw("COLZ"); - - cHit2->cd(5); - hXY5->Draw("COLZ"); - - cHit2->cd(6); - hXY6->Draw("COLZ"); - -} diff --git a/NPAnalysis/Gaspard/hit.c b/NPAnalysis/Gaspard/hit.c deleted file mode 100644 index 50db19da5e7c54b1cba977ffe34232fb328dcdd0..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/hit.c +++ /dev/null @@ -1,44 +0,0 @@ -{ - gROOT->Reset(); - gStyle->SetOptStat(1); - gStyle->SetPalette(51,0); - /////////////////////// - /////////////////////// -TFile *file0 = TFile::Open("./Result/myResult.root"); - - cEff = new TCanvas("cEff","Theta Distribution" ,100,100,600,600); - hTheta->Draw(); - - - cCM = new TCanvas("cCm" , "ThetaCM" , 100 , 100 , 600, 600) ; - hThetaCM->Draw(); - - cHit = new TCanvas("cHit","Hit" ,100,100,600,600); - hXY->Draw("COLZ"); - - - cHit2 = new TCanvas("cHit2","Hit" ,100,100,600,600); - cHit2->Divide(3,2); - - cHit2->cd(1); - hXY1->Draw("COLZ"); - - cHit2->cd(2); - hXY2->Draw("COLZ"); - - cHit2->cd(3); - hXY3->Draw("COLZ"); - - cHit2->cd(4); - hXY4->Draw("COLZ"); - - cHit2->cd(5); - hXY5->Draw("COLZ"); - - cHit2->cd(6); - hXY6->Draw("COLZ"); - - cEx = new TCanvas("cEx","Excitation Energy" ,100,100,300,300); - hEx->Draw(); - -} diff --git a/NPAnalysis/Gaspard/include/ConfigurationReader.hh b/NPAnalysis/Gaspard/include/ConfigurationReader.hh deleted file mode 100644 index d4b2e1fe1eaf30c39807391b226329724bab7b84..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/include/ConfigurationReader.hh +++ /dev/null @@ -1,20 +0,0 @@ -#ifndef ConfigurationReader_h -#define ConfigurationReader_h 1 - -// STL -#include <string> -#include <iostream> -#include <fstream> - -// NPL -#include "NPReaction.h" -using namespace NPL ; - -#include "must.hh" - -void ReadConfiguration ( string Path , array* myArray ) ; -Reaction* ReadReaction ( string Path ) ; - - -#endif - diff --git a/NPAnalysis/Gaspard/include/DetectorManager.hh b/NPAnalysis/Gaspard/include/DetectorManager.hh new file mode 100644 index 0000000000000000000000000000000000000000..eb49f01ef92509cd0531f4f12c7eace833047e61 --- /dev/null +++ b/NPAnalysis/Gaspard/include/DetectorManager.hh @@ -0,0 +1,39 @@ +#ifndef DetectorManager_h +#define DetectorManager_h + +// NPL +#include "VDetector.h" + +// STL +#include <string> +#include <map> + +using namespace std ; +using namespace NPA ; + +// This class manage a map of virtual detector + +class DetectorManager + { + public: + DetectorManager() ; + ~DetectorManager() ; + + public: + // Read stream at Path and pick-up Token declaration of Detector + void ReadConfigurationFile(string Path) ; + void BuildPhysicalEvent() ; + void BuildSimplePhysicalEvent() ; + void InitializeRootInput() ; + void InitializeRootOutput() ; + void AddDetector(string,VDetector*) ; + void ClearEventPhysics() ; + void ClearEventData() ; + + public: // The map containning all detectors + // Using a Map one can access to any detector using its name + map<string,VDetector*> m_Detector ; + + }; + +#endif diff --git a/NPAnalysis/Gaspard/include/ObjectManager.hh b/NPAnalysis/Gaspard/include/ObjectManager.hh index 9f90b2a5dc4c538d683f5bf378a2ece2250c591a..f7aa582ee7751c082d28923a2951e71ce38137f2 100644 --- a/NPAnalysis/Gaspard/include/ObjectManager.hh +++ b/NPAnalysis/Gaspard/include/ObjectManager.hh @@ -1,26 +1,65 @@ // You can use this file to declare your spectra, file, energy loss , ... and whatever you want. // This way you can remove all unnecessary declaration in the main programm. -// In order to help debugging and organising we use Name Space. +// In order to help debugging and organizing we use Name Space. + +///////////////////////////////////////////////////////////////////////////////////////////////// +// -------------------------------------- VARIOUS INCLUDE --------------------------------------- + +// NPA +#include "DetectorManager.hh" +#include "Must2Array.h" +#include "GaspardTracker.h" + +// STL C++ +#include <iostream> +#include <fstream> +#include <sstream> +#include <string> +#include <cmath> +#include <cstdlib> + +// ROOT +#include <TROOT.h> +#include <TChain.h> +#include <TFile.h> +#include <TLeaf.h> +#include <TVector3.h> +#include <TRandom.h> + +// NPL +#include "TMust2Data.h" +#include "TMust2Physics.h" +#include "NPReaction.h" +#include "RootInput.h" +#include "RootOutput.h" + +// Use CLHEP System of unit and Physical Constant +#include "CLHEP/Units/GlobalSystemOfUnits.h" +#include "CLHEP/Units/PhysicalConstants.h" -// -----------------------------------DOUBLE, INT, BOOL-------------------------------------------- - +// ---------------------------------------------------------------------------------------------- +double ThetaCalculation (TVector3 A , TVector3 B) ; +///////////////////////////////////////////////////////////////////////////////////////////////// +// ----------------------------------- DOUBLE, INT, BOOL AND MORE ------------------------------- namespace VARIABLE { - double Xs,Ys,Zs,N ; - double Theta = 0 ; - double Ex = 0 ; - double ThetaCM = 0 ; - int X,Y ; - bool check_light = false ; + // Declare your Variable here: + + double X1,Y1,Z1 ; + int N1,N2 = 0 ; + bool check= false ; + // A Usefull Simple Random Generator + TRandom Rand; } - + using namespace VARIABLE ; // ---------------------------------------------------------------------------------------------- +///////////////////////////////////////////////////////////////////////////////////////////////// // -----------------------------------GRAPH------------------------------------------------------ #include <TObject.h> #include <TH1.h> @@ -32,53 +71,48 @@ using namespace VARIABLE ; namespace GRAPH { // Declare your Spectra here: - const int NbElements = 18; - TH2F *hXY1 = new TH2F("hXY1","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY2 = new TH2F("hXY2","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY3 = new TH2F("hXY3","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY4 = new TH2F("hXY4","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY5 = new TH2F("hXY5","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY6 = new TH2F("hXY6","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hHitPattern[NbElements]; - - TH1F *hEnerX = new TH1F("hEnerX", "Energy X", 1000, -1, 10); - TH1F *hEnerY = new TH1F("hEnerY", "Energy Y", 1000, -1, 10); - TH2F *hEnerXY = new TH2F("hEnerXY", "Energy X vs Y", 100, -1, 10, 100, -1, 10); - TH2F *hEnerTheta = new TH2F("hEnerTheta", "Energy v.s. Theta", 90, 0, 180, 400, 0, 100); - - TH2F *hXY = new TH2F("hXY","", 300 , -130 , 130 , 300 , -130 , 130 ) ; - TH2F *hEA = new TH2F("hEA" ,"", 120 , 0 , 60 , 260 , 0 , 130 ) ; - TH1F *hEx = new TH1F("hEx", "E* (MeV)", 10000, -10, 40); - TH1F *hEHexa = new TH1F("hEHexa","EHexa (MeV)", 1000 , -5 , 5 ) ; - TH1F *hThetaCM = new TH1F("hThetaCM","ThetaCM (deg)", 30 , 0 , 60 ) ; - TH1F *hTheta = new TH1F("hTheta","Theta Lab (deg)",180,0,180) ; - TH2F *hThetaLabCM = new TH2F("hThetaLabCM","", 90 , 0 , 180 , 90 , 0 , 180 ) ; + TH1F *myHist1D = new TH1F("Hist1D","Histogramm 1D ; x ; count", 1000 , -5 , 5 ) ; + + TH2F *myHist2D = new TH2F("Hist2D","Histogramm 2D ; x ; y ", 128 , 1 , 128 , 128 , 1 , 128 ) ; + + } + +using namespace GRAPH ; +// -------------------------------------------------------------------------------------------- + + +/////////////////////////////////////////////////////////////////////////////////////////////// +// -----------------------------------CUT------------------------------------------------------ +#include <TCutG.h> +namespace CUT + { // Declare your Cut here: } -using namespace GRAPH ; -// ---------------------------------------------------------------------------------------------- +using namespace CUT ; +// -------------------------------------------------------------------------------------------- + -// -----------------------------------ENERGY LOSS------------------------------------------------ -/* -#include "NPEnergyLossCorrection.h" +//////////////////////////////////////////////////////////////////////////////////////////////// +// -----------------------------------ENERGY LOSS---------------------------------------------- +#include "NPEnergyLoss.h" using namespace NPL ; namespace ENERGYLOSS { -// EnergyLoss ProtonTarget = EnergyLoss ( "/home/Adrien/Desktop/dEdX.txt" , -// 100 , -// 1 ); - + // Declare your Energy loss here : + /* EnergyLoss ProtonTarget = EnergyLoss ( "CD2.txt" , + 100 , + 1 ); + */ } using namespace ENERGYLOSS ; -*/ // ---------------------------------------------------------------------------------------------- - +///////////////////////////////////////////////////////////////////////////////////////////////// diff --git a/NPAnalysis/Gaspard/include/must.hh b/NPAnalysis/Gaspard/include/must.hh deleted file mode 100644 index d3fa3a733fc029e7ed403f3a944e9dd68fc390c4..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/include/must.hh +++ /dev/null @@ -1,75 +0,0 @@ -#ifndef must_h -#define must_h 1 - -#include <vector> -#include <TFile.h> -#include <TH2.h> -#include <TH2F.h> -#include "CLHEP/Vector/ThreeVector.h" -using namespace std ; -using namespace CLHEP ; - -class must -{ -public: - must( double theta , double phi , double distance , double beta_u , double beta_v, double beta_w) ; - must( Hep3Vector C_X1_Y1 , - Hep3Vector C_X128_Y1 , - Hep3Vector C_X1_Y128 , - Hep3Vector C_X128_Y128 , - Hep3Vector TargetPosition ); - - double GetStripPositionX( int X , int Y ) ; - double GetStripPositionY( int X , int Y ) ; - double GetStripPositionZ( int X , int Y ) ; - void Print() ; - -private: - - vector< vector < double > > StripPositionX ; - vector< vector < double > > StripPositionY ; - vector< vector < double > > StripPositionZ ; - -}; - - -class array -{ -public: - array(); - void addTelescope( double theta , double phi , double distance , double beta_u , double beta_v , double beta_w) ; - void addTelescope( Hep3Vector C_X1_Y1 , - Hep3Vector C_X128_Y1 , - Hep3Vector C_X1_Y128 , - Hep3Vector C_X128_Y128 ); - - double GetStripPositionX( int N , int X , int Y ) ; - double GetStripPositionY( int N , int X , int Y ) ; - double GetStripPositionZ( int N , int X , int Y ) ; - double GetNumberOfTelescope() {return NumberOfTelescope;} ; - void Print() ; - -private: - vector< must > myArray ; - int NumberOfTelescope ; - -public: // Histogram - -/*// Hit Density - vector<TH2F*> HitDensity ; - - // Time of Flight Identification - vector<TH2F*> TOF ; - - // EDE Identification - vector<TH2F*> EDE ; - - // Kinematic Line - KinematicLine ; - - // Excitation Energy - ExcitationEnergy ;*/ - -}; - -#endif diff --git a/NPAnalysis/Gaspard/kinematic1MeV.contour b/NPAnalysis/Gaspard/kinematic1MeV.contour deleted file mode 100644 index 96c8267de29aec5472c1ec72fdbcf73171d62824..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/kinematic1MeV.contour +++ /dev/null @@ -1,27 +0,0 @@ -{ -//========= Macro generated from object: Graph/Graph -//========= by ROOT version5.22/00 - - TGraph *graph = new TGraph(9); - graph->SetName("Graph"); - graph->SetTitle("Graph"); - graph->SetFillColor(1); - graph->SetPoint(0,42.5558,25.636); - graph->SetPoint(1,36.4453,18.3991); - graph->SetPoint(2,21.2109,12.2588); - graph->SetPoint(3,6.22768,10.943); - graph->SetPoint(4,6.89732,5.67982); - graph->SetPoint(5,24.308,7.87281); - graph->SetPoint(6,38.8728,14.6711); - graph->SetPoint(7,44.6484,23.6623); - graph->SetPoint(8,42.6395,25.636); - - TH1 *Graph1 = new TH1F("Graph1","Graph",100,2.3856,48.4905); - Graph1->SetMinimum(3.68421); - Graph1->SetMaximum(27.6316); - Graph1->SetDirectory(0); - Graph1->SetStats(0); - graph->SetHistogram(Graph1); - - graph->Draw(""); -} diff --git a/NPAnalysis/Gaspard/makeAnalyse.sh b/NPAnalysis/Gaspard/makeAnalyse.sh deleted file mode 100755 index 6d62d9b5efef605ec9b14fafea319941f181dd9f..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/makeAnalyse.sh +++ /dev/null @@ -1,8 +0,0 @@ -# This file is here to help you making Analysis. -# Feel free to had more such as automatic lauch ou ROOT macro - - - -cd src ; make ; cd .. -./Analysis ../Inputs/Reaction/34Si.reaction ../Inputs/DetectorConfiguration/e530.detector -root -l hit.c diff --git a/NPAnalysis/Gaspard/root.ps b/NPAnalysis/Gaspard/root.ps deleted file mode 100644 index 79ba81bb0b06619d95cb69f587fb3b0ff93f94a4..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/root.ps +++ /dev/null @@ -1,98 +0,0 @@ -%!PS-Adobe-2.0 -%%Title: root.ps: figure d'article -%%Creator: ROOT Version 5.22/00 -%%CreationDate: Thu Sep 10 19:14:45 2009 -%%Orientation: Landscape -%%EndComments -%%BeginProlog -/s {stroke} def /l {lineto} def /m {moveto} def /t {translate} def -/sw {stringwidth} def /r {rotate} def /rl {roll} def /R {repeat} def -/d {rlineto} def /rm {rmoveto} def /gr {grestore} def /f {eofill} def -/c {setrgbcolor} def /black {0 setgray} def /sd {setdash} def -/cl {closepath} def /sf {scalefont setfont} def /lw {setlinewidth} def -/box {m dup 0 exch d exch 0 d 0 exch neg d cl} def -/NC{systemdict begin initclip end}def/C{NC box clip newpath}def -/bl {box s} def /bf {box f} def /Y { 0 exch d} def /X { 0 d} def -/mp {newpath /y exch def /x exch def} def -/side {[w .77 mul w .23 mul] .385 w mul sd w 0 l currentpoint t -144 r} def -/mr {mp x y w2 0 360 arc} def /m24 {mr s} def /m20 {mr f} def -/mb {mp x y w2 add m w2 neg 0 d 0 w neg d w 0 d 0 w d cl} def -/mt {mp x y w2 add m w2 neg w neg d w 0 d cl} def -/m21 {mb f} def /m25 {mb s} def /m22 {mt f} def /m26{mt s} def -/m23 {mp x y w2 sub m w2 w d w neg 0 d cl f} def -/m27 {mp x y w2 add m w3 neg w2 neg d w3 w2 neg d w3 w2 d cl s} def -/m28 {mp x w2 sub y w2 sub w3 add m w3 0 d 0 w3 neg d w3 0 d 0 w3 d w3 0 d 0 w3 d w3 neg 0 d 0 w3 d w3 neg 0 d - 0 w3 neg d w3 neg 0 d cl s } def -/m29 {mp gsave x w2 sub y w2 add w3 sub m currentpoint t - 4 {side} repeat cl fill gr} def -/m30 {mp gsave x w2 sub y w2 add w3 sub m currentpoint t - 4 {side} repeat cl s gr} def -/m31 {mp x y w2 sub m 0 w d x w2 sub y m w 0 d x w2 sub y w2 add m w w neg d x w2 sub y w2 - sub m w w d s} def -/m2 {mp x y w2 sub m 0 w d x w2 sub y m w 0 d s} def -/m5 {mp x w2 sub y w2 sub m w w d x w2 sub y w2 add m w w neg d s} def -/reEncode {exch findfont dup length dict begin {1 index /FID eq {pop pop} {def} ifelse } forall /Encoding exch def currentdict end dup /FontName get exch definefont pop } def [/Times-Bold /Times-Italic /Times-BoldItalic /Helvetica - /Helvetica-Oblique /Helvetica-Bold /Helvetica-BoldOblique /Courier /Courier-Oblique /Courier-Bold /Courier-BoldOblique /Times-Roman /AvantGarde-Book /AvantGarde-BookOblique /AvantGarde-Demi /AvantGarde-DemiOblique /Bookman-Demi - /Bookman-DemiItalic /Bookman-Light /Bookman-LightItalic /Helvetica-Narrow /Helvetica-Narrow-Bold /Helvetica-Narrow-BoldOblique /Helvetica-Narrow-Oblique /NewCenturySchlbk-Roman /NewCenturySchlbk-Bold - /NewCenturySchlbk-BoldItalic /NewCenturySchlbk-Italic /Palatino-Bold /Palatino-BoldItalic /Palatino-Italic /Palatino-Roman ] {ISOLatin1Encoding reEncode } forall -/oshow {gsave [] 0 sd true charpath stroke gr} def -/stwn { /fs exch def /fn exch def /text exch def fn findfont fs sf - text sw pop xs add /xs exch def} def -/stwb { /fs exch def /fn exch def /nbas exch def /textf exch deftextf length /tlen exch def nbas tlen gt {/nbas tlendef} iffn findfont fs sf textf dup length nbas sub nbas getinterval sw -pop neg xs add /xs exch def} def -/Zone {/iy exch def /ix exch def ix 1 sub 3144 mul 1 iy sub 2224 - mul t} def -%%EndProlog -%%BeginSetup -%%EndSetup -newpath gsave 90 r 0 -594 t 28 20 t .25 .25 scale gsave -%%Page: 1 1 - gsave gsave - 1 1 Zone - gsave 0 0 t black[ ] 0 sd 3 lw 0.999 0.999 0.999 c 2948 1992 0 0 bf black 0.999 0.999 0.999 c 2329 1574 324 219 bf black 2329 1574 324 219 bl 0.999 0.999 0.999 c 2329 1574 324 219 bf black 2329 1574 324 219 bl 1 1 1 c black 324 255 m 6 X 18 Y 5 X - -18 Y 6 X -18 Y 17 X -18 Y 5 X 36 Y 6 X -36 Y 6 X 36 Y 11 X -18 Y 11 X -18 Y 11 X 18 Y 17 X 36 Y 5 X -36 Y 23 X -18 Y 5 X 18 Y 6 X 36 Y 5 X -54 Y 6 X 36 Y 6 X -18 Y 5 X 18 Y 6 X -18 Y 17 X 18 Y 11 X -18 Y 11 X 53 Y 6 X -53 Y 5 X 18 Y 6 X -18 Y 11 X - 18 Y 11 X 18 Y 6 X 35 Y 5 X -18 Y 6 X -53 Y 5 X 18 Y 12 X 35 Y 5 X -53 Y 6 X -18 Y 5 X 107 Y 6 X -18 Y 5 X -53 Y 12 X 89 Y 5 X -89 Y 17 X -36 Y 6 X 89 Y 5 X -89 Y 6 X 71 Y 5 X -17 Y 6 X 17 Y 5 X 36 Y 6 X -18 Y 6 X 18 Y 5 X -53 Y 6 X 53 Y 5 X -36 Y - 6 X -53 Y 6 X 107 Y 5 X -89 Y 6 X -18 Y 5 X 36 Y 6 X -36 Y 5 X 71 Y 6 X -35 Y 6 X 71 Y 5 X 53 Y 6 X -53 Y 5 X -71 Y 6 X 53 Y 6 X -36 Y 5 X -17 Y 6 X -18 Y 5 X 125 Y 6 X -90 Y 6 X -17 Y 5 X 71 Y 6 X -89 Y 5 X 18 Y 6 X 17 Y 5 X 36 Y 6 X 36 Y 6 X -36 - Y 11 X 36 Y 5 X -125 Y 6 X 53 Y 6 X -35 Y 5 X 160 Y 6 X -160 Y 5 X 125 Y 6 X -72 Y 5 X -53 Y 6 X 35 Y 6 X 54 Y 5 X -54 Y 6 X 107 Y 5 X -124 Y 6 X 35 Y 6 X 54 Y 5 X -36 Y 6 X 18 Y 5 X 36 Y 6 X -90 Y 5 X -35 Y 6 X 53 Y 6 X -35 Y 5 X 35 Y 6 X 72 Y 5 X - -90 Y 6 X 54 Y 6 X -18 Y 5 X 18 Y 6 X -54 Y 5 X 36 Y 6 X 54 Y 6 X -107 Y 5 X 124 Y 6 X 18 Y 5 X -107 Y 6 X 54 Y 5 X -18 Y 6 X 36 Y 6 X 17 Y 11 X 72 Y 5 X -72 Y 6 X -107 Y 6 X 72 Y 5 X 53 Y 6 X 54 Y 5 X -89 Y 6 X -54 Y 5 X 71 Y 6 X -35 Y 6 X -54 Y 5 - X 18 Y 6 X 36 Y 5 X 143 Y 6 X 35 Y 6 X -89 Y 5 X -89 Y 6 X 53 Y 5 X 36 Y 6 X -54 Y 5 X 18 Y 6 X 143 Y 6 X -125 Y 5 X 18 Y 6 X 54 Y 5 X 53 Y 6 X 18 Y 6 X -89 Y 11 X 107 Y 5 X -36 Y 6 X -36 Y 5 X 72 Y 12 X 18 Y 5 X 71 Y 6 X -107 Y 5 X 214 Y 6 X -54 Y - 6 X -17 Y 5 X -36 Y 6 X -107 Y 5 X 232 Y 6 X -54 Y 6 X 107 Y 5 X -89 Y 6 X 214 Y 5 X -160 Y 6 X -72 Y 5 X 18 Y 6 X 232 Y 6 X 161 Y 5 X -179 Y 6 X 250 Y 5 X -71 Y 6 X -179 Y 6 X 143 Y 5 X 285 Y 6 X 36 Y 5 X -125 Y 6 X 357 Y 5 X -286 Y 6 X 215 Y 6 X - -286 Y 11 X 179 Y 5 X -143 Y 6 X -72 Y 6 X 250 Y 5 X -35 Y 6 X -72 Y 5 X 286 Y 6 X -90 Y 5 X -142 Y 12 X -304 Y 5 X 125 Y 6 X 107 Y 11 X -285 Y 6 X 196 Y 5 X -125 Y 6 X -53 Y 5 X -72 Y 6 X -160 Y 5 X 285 Y 6 X -410 Y 6 X 53 Y 5 X -124 Y 6 X 214 Y 5 - X -18 Y 6 X -36 Y 6 X -125 Y 5 X -160 Y 6 X 71 Y 5 X -18 Y 6 X 143 Y 6 X -143 Y 5 X -18 Y 6 X 36 Y 5 X -125 Y 6 X -18 Y 5 X 107 Y 6 X -71 Y 6 X 18 Y 5 X -18 Y 11 X 54 Y 6 X -107 Y 11 X -72 Y 6 X -36 Y 5 X 108 Y 6 X -36 Y 5 X 18 Y 6 X -18 Y 6 X 53 Y - 5 X -125 Y 11 X 54 Y 6 X 18 Y 6 X 18 Y 5 X -125 Y 6 X -18 Y 11 X -36 Y 5 X -18 Y 6 X 90 Y 11 X -18 Y 6 X 35 Y 5 X -17 Y 6 X -72 Y 6 X 54 Y 5 X 71 Y 6 X -53 Y 5 X -18 Y 6 X -18 Y 11 X -54 Y 6 X 54 Y 5 X -18 Y 6 X -36 Y 5 X -17 Y 12 X -18 Y 5 X 71 Y - 6 X -36 Y 5 X 18 Y 6 X 89 Y 6 X -107 Y 16 X 18 Y 6 X -53 Y 5 X 35 Y 12 X -53 Y 5 X 71 Y 6 X 18 Y 5 X -18 Y 6 X -18 Y 6 X 72 Y 5 X 18 Y 6 X -143 Y 5 X 89 Y 6 X -18 Y 5 X -53 Y 6 X 35 Y 6 X -17 Y 5 X -18 Y 6 X 35 Y 5 X -17 Y 6 X 53 Y 6 X -36 Y 5 X - -17 Y 6 X -54 Y 5 X 71 Y 6 X -17 Y 17 X -18 Y 5 X 35 Y 6 X 90 Y 5 X -125 Y 6 X 35 Y 6 X -35 Y 5 X 35 Y 6 X -17 Y 5 X -18 Y 6 X 53 Y 6 X -18 Y 5 X -35 Y 11 X 35 Y 6 X -53 Y 11 X 18 Y 6 X 35 Y 5 X -17 Y 6 X 17 Y 5 X 36 Y 17 X -53 Y 6 X -36 Y 5 X 53 Y - 6 X -17 Y 22 X 35 Y 6 X 54 Y 5 X -89 Y 12 X -18 Y 5 X 18 Y 6 X 53 Y 5 X -36 Y 6 X -17 Y 5 X 53 Y 17 X -89 Y 11 X 18 Y 12 X 89 Y 5 X -89 Y 6 X 18 Y 5 X -18 Y 6 X -18 Y 5 X 53 Y 6 X -53 Y 6 X 36 Y 5 X -36 Y 6 X 18 Y 5 X 18 Y 6 X -18 Y 6 X 53 Y 5 X - -18 Y 6 X -35 Y 11 X 35 Y 6 X -53 Y 5 X 53 Y 6 X -35 Y 11 X -18 Y 11 X -18 Y 6 X 54 Y 11 X -36 Y 5 X 36 Y 12 X -36 Y 11 X 36 Y 5 X -54 Y 6 X 71 Y 5 X -35 Y 6 X 18 Y 6 X 17 Y 5 X -53 Y 6 X 71 Y 5 X -71 Y 12 X 18 Y 5 X 18 Y 6 X -36 Y 5 X 89 Y 6 X -89 - Y 5 X 18 Y 6 X 18 Y 11 X s 324 219 m 2329 X s 369 266 m -47 Y s 481 243 m -24 Y s 592 243 m -24 Y s 704 243 m -24 Y s 816 243 m -24 Y s 927 266 m -47 Y s 1039 243 m -24 Y s 1151 243 m -24 Y s 1263 243 m -24 Y s 1374 243 m -24 Y s 1486 266 m -47 Y s - 1598 243 m -24 Y s 1709 243 m -24 Y s 1821 243 m -24 Y s 1933 243 m -24 Y s 2044 266 m -47 Y s 2156 243 m -24 Y s 2268 243 m -24 Y s 2380 243 m -24 Y s 2491 243 m -24 Y s 2603 266 m -47 Y s 369 266 m -47 Y s 2603 266 m -47 Y s - gsave 2948 1992 0 0 C 351 149 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (1) show NC gr - gsave 2948 1992 0 0 C 873 149 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (1.5) show NC gr - gsave 2948 1992 0 0 C 1461 149 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (2) show NC gr - gsave 2948 1992 0 0 C 1992 149 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (2.5) show NC gr - gsave 2948 1992 0 0 C 2580 149 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (3) show NC gr 324 1793 m 2329 X s 369 1745 m 48 Y s 481 1769 m 24 Y s 592 1769 m 24 Y s 704 1769 m 24 Y s 816 1769 m 24 Y s 927 1745 m 48 Y s 1039 1769 m 24 Y s 1151 - 1769 m 24 Y s 1263 1769 m 24 Y s 1374 1769 m 24 Y s 1486 1745 m 48 Y s 1598 1769 m 24 Y s 1709 1769 m 24 Y s 1821 1769 m 24 Y s 1933 1769 m 24 Y s 2044 1745 m 48 Y s 2156 1769 m 24 Y s 2268 1769 m 24 Y s 2380 1769 m 24 Y s 2491 1769 m 24 Y s 2603 - 1745 m 48 Y s 369 1745 m 48 Y s 2603 1745 m 48 Y s 324 219 m 1574 Y s 394 219 m -70 X s 359 255 m -35 X s 359 290 m -35 X s 359 326 m -35 X s 359 362 m -35 X s 394 397 m -70 X s 359 433 m -35 X s 359 469 m -35 X s 359 505 m -35 X s 359 540 m -35 X - s 394 576 m -70 X s 359 612 m -35 X s 359 647 m -35 X s 359 683 m -35 X s 359 719 m -35 X s 394 754 m -70 X s 359 790 m -35 X s 359 826 m -35 X s 359 861 m -35 X s 359 897 m -35 X s 394 933 m -70 X s 359 968 m -35 X s 359 1004 m -35 X s 359 1040 m - -35 X s 359 1075 m -35 X s 394 1111 m -70 X s 359 1147 m -35 X s 359 1182 m -35 X s 359 1218 m -35 X s 359 1254 m -35 X s 394 1289 m -70 X s 359 1325 m -35 X s 359 1361 m -35 X s 359 1396 m -35 X s 359 1432 m -35 X s 394 1468 m -70 X s 359 1504 m - -35 X s 359 1539 m -35 X s 359 1575 m -35 X s 359 1611 m -35 X s 394 1646 m -70 X s 394 1646 m -70 X s 359 1682 m -35 X s 359 1718 m -35 X s 359 1753 m -35 X s 359 1789 m -35 X s - gsave 2948 1992 0 0 C 268 189 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (0) show NC gr - gsave 2948 1992 0 0 C 224 369 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (10) show NC gr - gsave 2948 1992 0 0 C 224 548 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (20) show NC gr - gsave 2948 1992 0 0 C 224 724 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (30) show NC gr - gsave 2948 1992 0 0 C 224 904 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (40) show NC gr - gsave 2948 1992 0 0 C 224 1084 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (50) show NC gr - gsave 2948 1992 0 0 C 224 1259 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (60) show NC gr - gsave 2948 1992 0 0 C 224 1439 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (70) show NC gr - gsave 2948 1992 0 0 C 224 1619 t 0 r /Helvetica-Bold findfont 74.5782 sf 0 0 m (80) show NC gr 2653 219 m 1574 Y s 2583 219 m 70 X s 2618 255 m 35 X s 2618 290 m 35 X s 2618 326 m 35 X s 2618 362 m 35 X s 2583 397 m 70 X s 2618 433 m 35 X s 2618 - 469 m 35 X s 2618 505 m 35 X s 2618 540 m 35 X s 2583 576 m 70 X s 2618 612 m 35 X s 2618 647 m 35 X s 2618 683 m 35 X s 2618 719 m 35 X s 2583 754 m 70 X s 2618 790 m 35 X s 2618 826 m 35 X s 2618 861 m 35 X s 2618 897 m 35 X s 2583 933 m 70 X s - 2618 968 m 35 X s 2618 1004 m 35 X s 2618 1040 m 35 X s 2618 1075 m 35 X s 2583 1111 m 70 X s 2618 1147 m 35 X s 2618 1182 m 35 X s 2618 1218 m 35 X s 2618 1254 m 35 X s 2583 1289 m 70 X s 2618 1325 m 35 X s 2618 1361 m 35 X s 2618 1396 m 35 X s - 2618 1432 m 35 X s 2583 1468 m 70 X s 2618 1504 m 35 X s 2618 1539 m 35 X s 2618 1575 m 35 X s 2618 1611 m 35 X s 2583 1646 m 70 X s 2583 1646 m 70 X s 2618 1682 m 35 X s 2618 1718 m 35 X s 2618 1753 m 35 X s 2618 1789 m 35 X s 0.999 0.999 0.999 c - 432 111 29 1871 bf black 43 1871 m -9 Y 427 X 107 Y -9 X -98 Y f 29 1871 m 111 Y 432 X -111 Y -432 X cl s 1 1 1 c black - gsave 2948 1992 0 0 C 61 1904 t 0 r /Helvetica-Bold findfont 87.739 sf 0 0 m (E* \(MeV\)) show NC gr gr -showpage - gr -%%Trailer -%%Pages: 1 - gr gr gr -%%EOF diff --git a/NPAnalysis/Gaspard/src/Analysis.cc b/NPAnalysis/Gaspard/src/Analysis.cc index ad75bb183b1ccf0651bc89fa15d53c912c222f65..f3b8bbc888e5f02b4a3dcfe6a5e9610f8713a10a 100644 --- a/NPAnalysis/Gaspard/src/Analysis.cc +++ b/NPAnalysis/Gaspard/src/Analysis.cc @@ -1,337 +1,96 @@ -#include "must.hh" - -#include <fstream> -#include "sstream" -#include <string> -#include <cmath> - -//include ROOT -#include "TROOT.h" -#include "TChain.h" -#include "TFile.h" -#include "TString.h" -#include "TH2.h" -#include "TStyle.h" -#include "TRandom.h" - -//NPL -#include "TGaspardTrackerData.h" -#include "TInteractionCoordinates.h" -#include "TInitialConditions.h" -#include "NPReaction.h" -//#include "NPEnergyLossCorrection.h" -using namespace NPL ; - -// Use CLHEP System of unit and Physical Constant -#include "CLHEP/Units/GlobalSystemOfUnits.h" -#include "CLHEP/Units/PhysicalConstants.h" -#include "CLHEP/Vector/ThreeVector.h" - -//NPA -#include "ConfigurationReader.hh" #include "ObjectManager.hh" -// prototype -double ThetaCalculation( double StripX , double StripY ,double StripZ, double PPAC_X, double PPAC_Y, double PPAC_Theta); -double HexaEnergy(double HeavyEnergy ,double ThetaHeavy , double Excitation); -///////////// +using namespace std; -int main(int argc,char** argv) -{ - // A Usefull Random Generator - TRandom rand; - // Get arguments from command line +int main(int argc,char** argv) +{ + // test if number of arguments is correct if (argc != 4) { cout << - "you need to specify a Reaction file, a Detector file and a ROOT file such as : Analysis myReaction.reaction myDetector.detector myRootFile.root" - << endl; + "you need to specify both a Reaction file and a Detector file such as : Analysis myReaction.reaction myDetector.detector runToRead.run" + << endl; return 0; - } - string reactionfileName = argv[1]; - string detectorfileName = argv[2]; - string rootfileName = argv[3]; + } - ////////////////////////////////////////////////////////////////////// - ///////////////////////// Load Reaction ////////////////////////////// - ////////////////////////////////////////////////////////////////////// - Reaction* myReaction = ReadReaction(reactionfileName); + // get arguments + string reactionfileName = argv[1]; + string detectorfileName = argv[2]; + string runToReadfileName = argv[3]; - ////////////////////////////////////////////////////////////////////// - ///////////////////////// Load Detector ///////////////////////////// - ////////////////////////////////////////////////////////////////////// - array* myArray = new array(); - ReadConfiguration(detectorfileName, myArray); - - cout << "////////////////////////////////////////////" << endl; - cout << "Total Number Of telescope : " << myArray->GetNumberOfTelescope() << endl; - cout << "////////////////////////////////////////////" << endl; + // Instantiate RootInput and RootOutput singleton classes + RootInput:: getInstance(runToReadfileName); + RootOutput::getInstance("Analysis/Gaspard_AnalyzedData", "AnalyzedTree"); - ////////////////////////////////////////////////////////////////////// - /////////////////////// Load ROOT file /////////////////////////////// - ////////////////////////////////////////////////////////////////////// - // Open output ROOT file from NPTool simulation run - string path = getenv("NPTOOL"); - path += "/Outputs/Simulation/" + rootfileName; - TFile *inFile = new TFile(path.c_str()); - TTree *tree = (TTree*) inFile->Get("SimulatedTree"); + // Initialize the reaction + NPL::Reaction* myReaction = new Reaction(); + myReaction->ReadConfigurationFile(reactionfileName); -// TChain* t1 = new TChain("SimulatedTree"); -// t1->Add("../../Outputs/Simulation/mySimul.root"); - TGaspardTrackerData* EventGPD = new TGaspardTrackerData(); - tree->SetBranchAddress("GASPARD",&EventGPD); - TInteractionCoordinates* InterCoord = new TInteractionCoordinates(); - tree->SetBranchAddress("InteractionCoordinates",&InterCoord); - TInitialConditions* InitCond = new TInitialConditions(); - tree->SetBranchAddress("InitialConditions",&InitCond); + // Initialize the detector + DetectorManager* myDetector = new DetectorManager; + myDetector -> ReadConfigurationFile(detectorfileName); - ////////////////////////////////////////////////////////////////////// - /////////////////////// Load ROOT file /////////////////////////////// - ////////////////////////////////////////////////////////////////////// - for (Int_t i = 0; i < NbElements; i++) { - TString hname = Form("hHitPattern%d", i+1); - TString htitle = Form("HitPattern%d", i+1); - hHitPattern[i] = new TH2F(hname, htitle, 128, 1, 128, 128, 1, 128); - } + // Attach more branch to the output + double Ex = 0 ; double EE = 0 ; double TT = 0 ; double X = 0 ; double Y = 0 ; int det ; + RootOutput::getInstance()->GetTree()->Branch("ExcitationEnergy",&Ex,"Ex/D") ; + RootOutput::getInstance()->GetTree()->Branch("E",&EE,"EE/D") ; + RootOutput::getInstance()->GetTree()->Branch("A",&TT,"TT/D") ; + RootOutput::getInstance()->GetTree()->Branch("X",&X,"X/D") ; + RootOutput::getInstance()->GetTree()->Branch("Y",&Y,"Y/D") ; - //////////////// Analyse Part //////////////////////////////////////// - Int_t nentries = (Int_t) tree->GetEntries(); - cout << "There are " << nentries << " entries to process" << endl; - for (Int_t e = 0; e < nentries; e++) { - check_light =false ; - if ((e+1)%10000==0) cout << "Entries number " << e+1 << " on " << tree->GetEntries() << endl ; - tree->GetEntry(e); + // Get GaspardTracker pointer + GaspardTracker* GPDTrack = (GaspardTracker*) myDetector->m_Detector["GASPARD"]; - // first check that there is an event detected - if (EventGPD->GetGPDTrkFirstStageFrontEMult() > 0) { - // calculate multiplicity - Int_t multXE = EventGPD->GetGPDTrkFirstStageFrontEMult(); - Int_t multYE = EventGPD->GetGPDTrkFirstStageBackEMult(); - Int_t multXT = EventGPD->GetGPDTrkFirstStageFrontTMult(); - Int_t multYT = EventGPD->GetGPDTrkFirstStageBackTMult(); - // case of multiplicity 1 for the strips - if (multXE==1 && multYE==1 && multXT==1 && multYT==1) { - // calculate detector number - Int_t det_ref = EventGPD->GetGPDTrkFirstStageFrontEDetectorNbr(0); - Int_t detecXE = EventGPD->GetGPDTrkFirstStageFrontEDetectorNbr(0) / det_ref; - Int_t detecXT = EventGPD->GetGPDTrkFirstStageFrontTDetectorNbr(0) / det_ref; - Int_t detecYE = EventGPD->GetGPDTrkFirstStageBackEDetectorNbr(0) / det_ref; - Int_t detecYT = EventGPD->GetGPDTrkFirstStageBackTDetectorNbr(0) / det_ref; - det_ref -= 1000; // for TGaspardTrackerDummyShape - // case of same detector - if (detecXE*detecXT*detecYE*detecYT == 1) { - // calculate strip number - Int_t stripXE = EventGPD->GetGPDTrkFirstStageFrontEStripNbr(0); - Int_t stripXT = EventGPD->GetGPDTrkFirstStageFrontTStripNbr(0); - Int_t stripYE = EventGPD->GetGPDTrkFirstStageBackEStripNbr(0); - Int_t stripYT = EventGPD->GetGPDTrkFirstStageBackTStripNbr(0); - // case of same strips on X and Y - if (stripXE == stripXT && stripYE == stripYT) { // here we have a good strip event - // fill hit pattern histograms - hHitPattern[det_ref-1]->Fill(stripXE, stripYE); - // get energy from strips - Double_t EnergyStripX = EventGPD->GetGPDTrkFirstStageFrontEEnergy(0) * MeV; - Double_t EnergyStripY = EventGPD->GetGPDTrkFirstStageBackEEnergy(0) * MeV; - Double_t EnergyStrip = 0.5 * (EnergyStripX + EnergyStripY); -// Double_t EnergyStrip = EnergyStripX; -// if (EnergyStripY > EnergyStrip) EnergyStrip = EnergyStripY; - Double_t TotalEnergy = EnergyStrip; - // fill detected energy - hEnerX->Fill(EnergyStripX); - hEnerY->Fill(EnergyStripY); - hEnerXY->Fill(EnergyStripX, EnergyStripY); - // calculate multiplicity of 2nd and third stages - Int_t mult2E = EventGPD->GetGPDTrkSecondStageEMult(); - Int_t mult2T = EventGPD->GetGPDTrkSecondStageTMult(); - Int_t mult3E = EventGPD->GetGPDTrkThirdStageEMult(); - Int_t mult3T = EventGPD->GetGPDTrkThirdStageTMult(); - // check if we have a 2nd stage event - if (mult2E==1 && mult2T==1) { - Double_t EnergySecond = EventGPD->GetGPDTrkSecondStageEEnergy(0); - TotalEnergy += EnergySecond; - } - // check if we have a third stage event - if (mult3E==1 && mult3T==1) { - Double_t EnergyThird = EventGPD->GetGPDTrkThirdStageEEnergy(0); - TotalEnergy += EnergyThird; - } - // get emitted angle in the beam frame -// Double_t Theta = InitCond->GetICEmittedAngleThetaLabIncidentFrame(0) * deg; - Double_t Theta = InterCoord->GetDetectedAngleTheta(0) * deg; - hTheta->Fill(Theta / deg); - hEnerTheta->Fill(Theta / deg, TotalEnergy); - // calculate lab angle from strips - Double_t xstrip = myArray->GetStripPositionX(det_ref, stripXE, stripYE); - Double_t ystrip = myArray->GetStripPositionY(det_ref, stripXE, stripYE); - Double_t zstrip = myArray->GetStripPositionZ(det_ref, stripXE, stripYE); - Hep3Vector posstrip(xstrip, ystrip, zstrip); - Double_t ThetaStrip = posstrip.theta() * rad; - cout << Theta/deg << " " << ThetaStrip/deg << endl; - // calculate excitation energy -// Ex = myReaction->ReconstructRelativistic(TotalEnergy / MeV, Theta / rad); - Ex = myReaction->ReconstructRelativistic(TotalEnergy / MeV, ThetaStrip / rad); - // condition on c.m. angle - // for (d,p) reaction theta_c.m. < 40 deg - Double_t ThetaCM = myReaction->EnergyLabToThetaCM(TotalEnergy / MeV, -500) * rad; - hThetaCM->Fill(ThetaCM / deg); - hThetaLabCM -> Fill(ThetaCM/deg,Theta/deg); -// if (ThetaCM/deg < 20) { - hEx->Fill(Ex); -// } - } - else { - cout << "Not same strips" << endl; - } - } - else { - cout << "Not same detector" << endl; - } - } - else { - cout << "Multiplicity is not one, it is: " << endl; - cout << "\tmultXE: " << multXE << endl; - cout << "\tmultXT: " << multXT << endl; - cout << "\tmultYE: " << multYE << endl; - cout << "\tmultYT: " << multYT << endl; - } - } -/* - PhysicsMUST2->BuildPhysicalEvent(EventMUST2) ; - - /////////////////////// - // Put your Analysis Here - - if( PhysicsMUST2->TelescopeNumber.size()>0) - { - for(Int_t mult = 0 ; mult < PhysicsMUST2->TelescopeNumber.size() ; mult++) - { - N=PhysicsMUST2->TelescopeNumber.at(mult); - - Double_t Theta = rad * ThetaCalculation( - Xs=myArray->GetStripPositionX( N , X=PhysicsMUST2->Si_X.at(mult) , Y=PhysicsMUST2->Si_Y.at(mult) ) , - Ys=myArray->GetStripPositionY( N , X=PhysicsMUST2->Si_X.at(mult) , Y=PhysicsMUST2->Si_Y.at(mult) ) , - Zs=myArray->GetStripPositionZ( N , X=PhysicsMUST2->Si_X.at(mult) , Y=PhysicsMUST2->Si_Y.at(mult) ) , - BeamX , - BeamY , - BeamTheta ); - - // Rand prevent binning effect due to discret Hit position pattern (strip center) - hXY->Fill( (Xs-1)+2*rand.Uniform() ,(Ys-1)+2*rand.Uniform() ); - - if(N==1)hXY1 -> Fill(X,Y) ; - if(N==2)hXY2 -> Fill(X,Y) ; - if(N==3)hXY3 -> Fill(X,Y) ; - if(N==4)hXY4 -> Fill(X,Y) ; - - Double_t Energy = PhysicsMUST2->TotalEnergy.at(mult)*MeV ; - - // Correct Energy Loss in the target -// Energy = ProtonTarget.EvaluateInitialEnergy( Energy , -// 20.618556701*micrometer , -// Theta ); - - - Ex = myReaction->ReconstructRelativistic( Energy , Theta ) ; - - ThetaCM = myReaction->EnergyLabToThetaCM(Energy) ; - hThetaLabCM -> Fill(ThetaCM/deg,Theta/deg) ; - hThetaCM -> Fill(ThetaCM/deg) ; - hEA-> Fill(Theta/deg , Energy) ; - //hEA-> Fill(InitialTheta , InitialEnergy) ; - hTheta->Fill(Theta/deg) ; - check_light = true ; - - hEx->Fill(Ex); -*/ - /*if(N==5 && check_light) - { - - double Hexa = HexaEnergy(PhysicsMUST2->TotalEnergy.at(mult) ,Ex) ; - cout << Hexa << endl ; - hEHexa->Fill(Hexa); - }*/ -/* } - - } - - PhysicsMUST2->Clear(); - ////////////////////// - - PhysicsMUST2->Clear();*/ - } - - TFile outfile("./Result/myResult.root","RECREATE") ; - - for (Int_t i = 0; i < NbElements; i++) hHitPattern[i]->Write(); - hTheta -> Write() ; - hThetaCM -> Write() ; - hThetaLabCM -> Write() ; - hXY -> Write() ; - hXY1 -> Write() ; - hXY2 -> Write() ; - hXY3 -> Write() ; - hXY4 -> Write() ; - hXY5 -> Write() ; - hXY6 -> Write() ; - hEx -> Write() ; - hEA -> Write() ; - hEHexa -> Write() ; - hEnerX->Write(); - hEnerY->Write(); - hEnerXY->Write(); - hEnerTheta->Write(); + // Get the TChain and treat it + TChain* chain = RootInput:: getInstance() -> GetChain(); - outfile.Close() ; - -return 0 ; -} + // Analysis is here! + int nentries = chain->GetEntries(); + cout << "Number of entries to be analysed: " << nentries << endl; + + for (int i = 0; i < nentries; i ++) { + if (i%10000 == 0 && i!=0) cout << i << " analysed events" << endl; + chain -> GetEntry(i); + + // Treat Gaspard event + myDetector->ClearEventPhysics(); + myDetector->BuildPhysicalEvent(); + // Get total energy and coordinates of interaction + double E = GPDTrack->GetEnergyDeposit(); + TVector3 A = GPDTrack->GetPositionOfInteraction(); + // Calculate scattering angle + double Theta = ThetaCalculation (A ,TVector3(0,0,1)); -double ThetaCalculation( double StripX , double StripY ,double StripZ, double PPAC_X, double PPAC_Y, double PPAC_Theta) - { - double TargetZ = -50;//mm - StripX = StripX-PPAC_X ; StripY = StripY-PPAC_Y ; StripZ=StripZ-TargetZ ; - double norme = sqrt ( StripX*StripX + StripY*StripY + StripZ*StripZ ) ; - double Theta = acos ( StripZ / norme ) ; - return (Theta-PPAC_Theta) ; //return Theta in rad - } + // Calculate excitatioin energy + if (E > -1000) Ex = myReaction->ReconstructRelativistic(E, Theta); + else Ex = -100; + EE = E ; TT = Theta/deg; + if (E>-1000) { + X = A . X(); + Y = A . Y(); + } + else { + X = -1000 ; Y = -1000; + } -double HexaEnergy(double HeavyEnergy ,double ThetaHeavy , double Excitation) - { - /*double M10He = 9362.73 + Excitation ; - double Malpha = 3727.42 ; - double Mn = 939.56536 ; - double M2 = M10He * M10He ; - double Ma2 = Malpha * Malpha ; + // Fill output tree + RootOutput::getInstance()->GetTree()->Fill(); + } - // HeavyEnergyLab to HeavyEnergyCM - - // Initial Lab to CM Momentum - G4double InitialMomentum = sqrt(EnergyHeavy*EnergyHeavy - Ma2) ; - G4double InitialMomentumX = sin(ThetaHeavy) ; - G4double InitialMomentumY = sin(ThetaHeavy) ; - G4double InitialMomentumZ = cos(ThetaHeavy) ; - - // Beta and Gamma CM to Lab - G4double betaX = (DaughterMomentumX*p - InitialMomentumX*InitialMomentum)/E ; - G4double betaY = (DaughterMomentumY*p - InitialMomentumY*InitialMomentum)/E ; - G4double betaZ = (DaughterMomentumZ*p - InitialMomentumZ*InitialMomentum)/E ; - G4double beta = sqrt (betaX*betaX + betaY*betaY + betaZ*betaZ ) ; - G4double beta2 = beta*beta ; - G4double gamma = 1 / ( sqrt(1 - beta2 ) ) ; - // Calculate everything for heavy particule - - G4double NewEnergy = - gamma*E - - betaX*gamma*DaughterMomentumX*p - - betaY*gamma*DaughterMomentumY*p - - betaZ*gamma*DaughterMomentumZ*p; - - double HexaEnergy = sqrt( ( (M2 + Ma2)/(2*M10He) - (HeavyEnergy) ) * (2*M10He) ); - - return(6*Mn-HexaEnergy) ;*/ - } + // delete singleton classes + RootOutput::getInstance()->Destroy(); + RootInput::getInstance()->Destroy(); + return 0; +} + + +double ThetaCalculation (TVector3 A , TVector3 B) +{ + double Theta = acos( (A.Dot(B)) / (A.Mag()*B.Mag()) ); + return Theta ; +} diff --git a/NPAnalysis/Gaspard/src/ConfigurationReader.cc b/NPAnalysis/Gaspard/src/ConfigurationReader.cc deleted file mode 100644 index 827855f119ecf7c4073a5b1e828cfcf7d4fe3184..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/src/ConfigurationReader.cc +++ /dev/null @@ -1,319 +0,0 @@ -#include "ConfigurationReader.hh" -#include <vector> -#include <stdlib.h> -#include "CLHEP/Vector/ThreeVector.h" - -using namespace std ; -using namespace CLHEP ; -//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... -void ReadConfiguration( string Path , array* myArray) -{ - -////////General Reading needs//////// - string LineBuffer; - string DataBuffer; - -////////MUST2 Detector needs//////// - Hep3Vector A,B,C,D; - - double Theta,Phi,Distance,Beta_u,Beta_v,Beta_w; - int SI,SILI,CSI,VIS; - bool check_A=false;bool check_B=false;bool check_C=false;bool check_D=false; - bool check_Theta=false;bool check_Phi=false;bool check_Distance=false;bool check_Beta=false; - - bool MUST2=false; - bool Target=false; - bool CryoTarget=false; - -//////////Target Definition needs/////// - double TargetThickness, TargetRadius; - - //for cryogenic target - - double WindowsThickness, TargetTemperature, TargetPressure; - -////////////////////////////////////////////////////////////////////////////////////////// - ifstream ConfigFile; - ConfigFile.open(Path.c_str()); - - if( ConfigFile.is_open() ) - cout << " Configuration file " << Path << " loading " << endl; - else{ - cout << " Error, no configuration file" << Path << " found" << endl; return;} - - int i=0; - while( !ConfigFile.eof() ) - { - - //Pick-up next line - getline(ConfigFile, LineBuffer); i++; - - //Search for comment Symbol % - if(LineBuffer.compare(0,1,"%")==0) {;} - - //Search for Telescope - else if(LineBuffer.compare(0,9,"GPDSquare")==0 || LineBuffer.compare(0,13,"GPDDummyShape")==0) - { - - //A MUST2 Telescope is detected: - - //Check if it is the first telescope added. If yes, Instantiate a new Array linked to the logical World. - if(!MUST2) - { - //MUST2Array= new MUST2Array(world); - MUST2=true; - } - - check_A=false; check_B=false; check_C=false; check_D=false; - check_Theta=false; check_Phi=false; check_Distance=false; - cout << "///////////////" << endl; - cout << "Telescope found" << endl; - cout << "///////////////" << endl; - - ConfigFile >> DataBuffer; - //Position method - if(DataBuffer.compare(0,6,"X1_Y1=")==0) - { - check_A=true; - ConfigFile >> DataBuffer ; A.setX( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; A.setY( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; A.setZ( atof( DataBuffer.c_str() ) ); - - cout << "X1 Y1 corner position : " << A << endl; - } - - //Angle method - if(DataBuffer.compare(0,6,"THETA=")==0) - { - check_Theta=true; - ConfigFile >> DataBuffer ; Theta=atof(DataBuffer.c_str()) ; - cout << "Theta: "<< Theta << endl; - } - - - ConfigFile >> DataBuffer; - - //Position method - if(DataBuffer.compare(0,8,"X128_Y1=")==0) - { - check_B=true; - ConfigFile >> DataBuffer ; B.setX( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; B.setY( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; B.setZ( atof( DataBuffer.c_str() ) ); - - cout << "X128 Y1 corner position : " << B << endl; - } - - //Angle method - if(DataBuffer.compare(0,4,"PHI=")==0) - { - check_Phi=true; - ConfigFile >> DataBuffer ; Phi=atof(DataBuffer.c_str()) ; - cout << "Phi: "<< Phi << endl; - } - - ConfigFile >> DataBuffer; - - //Position method - if(DataBuffer.compare(0,8,"X1_Y128=")==0) - { - check_C=true; - ConfigFile >> DataBuffer ; C.setX( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; C.setY( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; C.setZ( atof( DataBuffer.c_str() ) ); - - cout << "X1 Y128 corner position : " << C << endl; - } - - //Angle method - if(DataBuffer.compare(0,2,"R=")==0) - { - check_Distance=true; - ConfigFile >> DataBuffer ; Distance=atof(DataBuffer.c_str()) ; - cout << "Distance: "<< Distance << endl; - } - - ConfigFile >> DataBuffer; - //Position method - if(DataBuffer.compare(0,10,"X128_Y128=")==0) - { - check_D=true; - ConfigFile >> DataBuffer ; D.setX( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; D.setY( atof( DataBuffer.c_str() ) ); - ConfigFile >> DataBuffer ; D.setZ( atof( DataBuffer.c_str() ) ); - - cout << "X128 Y128 corner position : " << D << endl; - } - - //Angle method - if(DataBuffer.compare(0,5,"BETA=")==0) - { - check_Beta=true; - ConfigFile >> DataBuffer ; Beta_u=atof(DataBuffer.c_str()) ; - ConfigFile >> DataBuffer ; Beta_v=atof(DataBuffer.c_str()) ; - ConfigFile >> DataBuffer ; Beta_w=atof(DataBuffer.c_str()) ; - cout << "Beta: "<< Beta_u << Beta_v << Beta_w << endl; - } - - ConfigFile >> DataBuffer; - if(DataBuffer.compare(0,10,"FIRSTSTAGE=")==0) - { - ConfigFile >> DataBuffer;SI=atof(DataBuffer.c_str()) ; - } - - ConfigFile >> DataBuffer; - if(DataBuffer.compare(0,11,"SECONDSTAGE=")==0) - { - ConfigFile >> DataBuffer;SILI=atof(DataBuffer.c_str()) ; - } - - - ConfigFile >> DataBuffer; - if(DataBuffer.compare(0,10,"THIRDSTAGE=")==0) - { - ConfigFile >> DataBuffer;CSI=atof(DataBuffer.c_str()) ; - } - - - ///Add The previously define telescope - - //With position method - if(check_A && check_B && check_C && check_D) - { - myArray->addTelescope( A , - B , - C , - D ); - - } - - //with angle method - if(check_Theta && check_Phi && check_Distance && check_Beta) - { - myArray->addTelescope( Theta , - Phi , - Distance , - Beta_u , - Beta_v , - Beta_w ); - - } - - - } - - //Nothing understandable - else ; - - - } - -ConfigFile.close(); -return; -} - -//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... -Reaction* ReadReaction( string Path ) -{ -////////General Reading needs//////// - string LineBuffer; - string DataBuffer; - -////////Reaction Setting needs/////// -string Beam,Target,Heavy,Light,CrossSectionPath; -double BeamEnergy,ExcitationEnergy; -Reaction* myReaction; - -////////////////////////////////////////////////////////////////////////////////////////// - ifstream ReactionFile; - ReactionFile.open(Path.c_str()); - - if( ReactionFile.is_open() ) - cout << " Reaction file " << Path << " loading " << endl; - else{ - cout << " Error, no Reaction file" << Path << " found" << endl; return myReaction;} - - - int i=0; - while( !ReactionFile.eof() ) - { - //Pick-up next line - getline(ReactionFile, LineBuffer); i++; - - //Search for comment Symbol % - if(LineBuffer.compare(0,1,"%")==0) {;} - - else if(LineBuffer.compare(0,9,"Transfert")==0) - { - - ReactionFile >> DataBuffer; - if(DataBuffer.compare(0,5,"Beam=")==0) - { - ReactionFile >> DataBuffer; - Beam = DataBuffer; - cout << "Beam " << Beam << endl; - } - - ReactionFile >> DataBuffer; - if(DataBuffer.compare(0,7,"Target=")==0) - { - ReactionFile >> DataBuffer; - Target = DataBuffer; - cout << "Target " << Target << endl; - } - - ReactionFile >> DataBuffer; - if(DataBuffer.compare(0,6,"Light=")==0) - { - ReactionFile >> DataBuffer; - Light = DataBuffer; - cout << "Light " << Light << endl; - } - - ReactionFile >> DataBuffer; - if(DataBuffer.compare(0,6,"Heavy=")==0) - { - ReactionFile >> DataBuffer; - Heavy = DataBuffer; - cout << "Heavy " << Heavy << endl; - } - - ReactionFile >> DataBuffer; - if(DataBuffer.compare(0,17,"ExcitationEnergy=")==0) - { - ReactionFile >> DataBuffer; - ExcitationEnergy = atof(DataBuffer.c_str()); - cout << "ExcitationEnergy " << ExcitationEnergy << " MeV" << endl; - } - - ReactionFile >> DataBuffer; - if(DataBuffer.compare(0,11,"BeamEnergy=")==0) - { - ReactionFile >> DataBuffer; - BeamEnergy = atof(DataBuffer.c_str()); - cout << "EnergyBeam " << BeamEnergy << " MeV" << endl; - } - - ReactionFile >> DataBuffer;ReactionFile >> DataBuffer; - ReactionFile >> DataBuffer;ReactionFile >> DataBuffer; - ReactionFile >> DataBuffer;ReactionFile >> DataBuffer; - ReactionFile >> DataBuffer;ReactionFile >> DataBuffer; - ReactionFile >> DataBuffer;ReactionFile >> DataBuffer; - - - ReactionFile >> DataBuffer; - if(DataBuffer.compare(0,17,"CrossSectionPath=")==0) - { - ReactionFile >> CrossSectionPath; - cout << "CrossSectionPath " << CrossSectionPath << endl; - } - } - } - -myReaction = new Reaction(Beam, Target, Light, Heavy,BeamEnergy,ExcitationEnergy,CrossSectionPath); - -ReactionFile.close(); -return myReaction; - -} - diff --git a/NPAnalysis/Gaspard/src/ConfigurationReader.hh b/NPAnalysis/Gaspard/src/ConfigurationReader.hh deleted file mode 100644 index d4b2e1fe1eaf30c39807391b226329724bab7b84..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/src/ConfigurationReader.hh +++ /dev/null @@ -1,20 +0,0 @@ -#ifndef ConfigurationReader_h -#define ConfigurationReader_h 1 - -// STL -#include <string> -#include <iostream> -#include <fstream> - -// NPL -#include "NPReaction.h" -using namespace NPL ; - -#include "must.hh" - -void ReadConfiguration ( string Path , array* myArray ) ; -Reaction* ReadReaction ( string Path ) ; - - -#endif - diff --git a/NPAnalysis/Gaspard/src/DetectorManager.cc b/NPAnalysis/Gaspard/src/DetectorManager.cc new file mode 100644 index 0000000000000000000000000000000000000000..93bbbd699e84eaacedecc00f80499a4608e9a676 --- /dev/null +++ b/NPAnalysis/Gaspard/src/DetectorManager.cc @@ -0,0 +1,229 @@ +#include "DetectorManager.hh" + +// STL +#include <iostream> +#include <fstream> +#include <cstdlib> + +// Detector +#include "Must2Array.h" +#include "GaspardTracker.h" + +///////////////////////////////////////////////////////////////////////////////////////////////// +// Default Constructor and Destructor +DetectorManager::DetectorManager() + {} + +///////////////////////////////////////////////////////////////////////////////////////////////// +DetectorManager::~DetectorManager() + {} + + +///////////////////////////////////////////////////////////////////////////////////////////////// +// Read stream at ConfigFile and pick-up Token declaration of Detector +void DetectorManager::ReadConfigurationFile(string Path) + { + ////////General Reading needs//////// + string LineBuffer; + string DataBuffer; + + /////////Boolean//////////////////// + bool MUST2 = false ; + bool AddThinSi = false ; + bool GeneralTarget = false ; + bool GPDTracker = false ; // Gaspard Tracker + ////////////////////////////////////////////////////////////////////////////////////////// + // added by Nicolas [07/05/09] + string GlobalPath = getenv("NPTOOL"); + Path = GlobalPath + "/Inputs/DetectorConfiguration/" + Path; + ifstream ConfigFile; + ConfigFile.open(Path.c_str()); + + if (ConfigFile.is_open()) + { + cout << "/////////////////////////////" << endl; + cout << " Configuration file " << Path << " loading " << endl; + } + + else + { + cout << " Error, no configuration file" << Path << " found" << endl; + return; + } + + + while (!ConfigFile.eof()) { + //Pick-up next line + getline(ConfigFile, LineBuffer); + //Search for comment Symbol: % + if (LineBuffer.compare(0, 1, "%") == 0) { /*Do Nothing*/ ;} + + //////////////////////////////////////////// + //////////// Search for Gaspard //////////// + //////////////////////////////////////////// + else if (LineBuffer.compare(0, 14, "GaspardTracker") == 0 && GPDTracker == false) { + GPDTracker = true ; + cout << "//////// Gaspard Tracker ////////" << endl ; + + // Instantiate the new array as a VDetector Object + VDetector* myDetector = new GaspardTracker() ; + + // Read Position of Telescope + ConfigFile.close() ; + myDetector->ReadConfiguration(Path) ; + ConfigFile.open(Path.c_str()) ; + + // Add array to the VDetector Vector + AddDetector("GASPARD", myDetector) ; + } + + //////////////////////////////////////////// + //////// Search for MUST2 Array //////// + //////////////////////////////////////////// + else if (LineBuffer.compare(0, 10, "MUST2Array") == 0 && MUST2 == false) { + MUST2 = true ; + cout << "//////// MUST2 Array ////////" << endl << endl ; + + // Instantiate the new array as a VDetector Object + VDetector* myDetector = new MUST2Array() ; + + // Read Position of Telescope + ConfigFile.close() ; + myDetector->ReadConfiguration(Path) ; + ConfigFile.open(Path.c_str()) ; + + // Add array to the VDetector Vector + AddDetector("MUST2" , myDetector) ; + } + + /* //////////////////////////////////////////// + ////////// Search for Add.ThinSi /////////// + //////////////////////////////////////////// + else if (LineBuffer.compare(0, 9, "AddThinSi") == 0 && AddThinSi == false) { + AddThinSi = true ; + cout << "//////// Thin Si ////////" << endl << endl ; + + // Instantiate the new array as a VDetector Object + VDetector* myDetector = new ThinSi() ; + + // Read Position of Telescope + ConfigFile.close() ; + myDetector->ReadConfiguration(Path) ; + ConfigFile.open(Path.c_str()) ; + + // Add array to the VDetector Vector + AddDetector(myDetector) ; + } + + //////////////////////////////////////////// + //////////// Search for Target ///////////// + //////////////////////////////////////////// + + else if (LineBuffer.compare(0, 13, "GeneralTarget") == 0 && GeneralTarget == false) { + GeneralTarget = true ; + cout << "////////// Target ///////////" << endl << endl ; + + // Instantiate the new array as a VDetector Objects + VDetector* myDetector = new Target() ; + + // Read Position and target specification + ConfigFile.close() ; + myDetector->ReadConfiguration(Path) ; + ConfigFile.open(Path.c_str()) ; + + m_TargetThickness = ((Target*)myDetector)->GetTargetThickness() ; + m_TargetRadius = ((Target*)myDetector)->GetTargetRadius() ; + m_TargetX = ((Target*)myDetector)->GetTargetX() ; + m_TargetY = ((Target*)myDetector)->GetTargetY() ; + m_TargetZ = ((Target*)myDetector)->GetTargetZ() ; + + // Add target to the VDetector Vector + AddDetector(myDetector) ; + }*/ + + //Nothing understandable + //else ; + } + + ConfigFile.close(); + + InitializeRootInput(); + InitializeRootOutput(); + + return ; + } + +///////////////////////////////////////////////////////////////////////////////////////////////// + +void DetectorManager::BuildPhysicalEvent() + { + map<string,VDetector*>::iterator it ; + + for( it = m_Detector.begin() ; it != m_Detector.end() ; ++it) + { + it->second->BuildPhysicalEvent() ; + } + } +///////////////////////////////////////////////////////////////////////////////////////////////// + +void DetectorManager::BuildSimplePhysicalEvent() + { + map<string,VDetector*>::iterator it ; + + for( it = m_Detector.begin() ; it != m_Detector.end() ; ++it) + { + it->second->BuildSimplePhysicalEvent() ; + } + } +///////////////////////////////////////////////////////////////////////////////////////////////// + +void DetectorManager::InitializeRootInput() + { + map<string,VDetector*>::iterator it ; + + for( it = m_Detector.begin() ; it != m_Detector.end() ; ++it) + { + it->second->InitializeRootInput() ; + } + } +///////////////////////////////////////////////////////////////////////////////////////////////// + +void DetectorManager::InitializeRootOutput() + { + map<string,VDetector*>::iterator it ; + + for( it = m_Detector.begin() ; it != m_Detector.end() ; ++it) + { + it->second->InitializeRootOutput() ; + } + + } +///////////////////////////////////////////////////////////////////////////////////////////////// + +void DetectorManager::AddDetector(string DetectorName , VDetector* newDetector) + { + m_Detector[DetectorName] = newDetector ; + } +///////////////////////////////////////////////////////////////////////////////////////////////// +void DetectorManager::ClearEventPhysics() + { + map<string,VDetector*>::iterator it ; + + for( it = m_Detector.begin() ; it != m_Detector.end() ; ++it) + { + it->second->ClearEventPhysics() ; + } + + } +///////////////////////////////////////////////////////////////////////////////////////////////// +void DetectorManager::ClearEventData() + { + map<string,VDetector*>::iterator it ; + + for( it = m_Detector.begin() ; it != m_Detector.end() ; ++it) + { + it->second->ClearEventData() ; + } + + } + diff --git a/NPAnalysis/Gaspard/src/GNUmakefile b/NPAnalysis/Gaspard/src/GNUmakefile index 51c5fef9d17ce76ec2688e1b0f19b11c285f5ab8..c18b108aaeff1d074e76e615e2f38d050e4cfedc 100644 --- a/NPAnalysis/Gaspard/src/GNUmakefile +++ b/NPAnalysis/Gaspard/src/GNUmakefile @@ -5,7 +5,7 @@ CPP=g++ EXEC=Analysis # local includes -NPAINCLUDES = $(NPANA)/include +NPAINCLUDES = ../include # ROOT includes CXXFLAGS += `root-config --cflags` @@ -15,21 +15,21 @@ CXXFLAGS += -I$(CLHEP_INCLUDE_DIR) CXXFLAGS += -I$(NPAINCLUDES) CXXFLAGS += -I$(NPLIB)/include - LDFLAGS = `root-config --libs` -lMathMore -#LDFLAGS += -L$(NPLIB) -lMust2Data -lReaction -lEnergyLossCorrection -LDFLAGS += -L$(NPLIB)/lib -lMust2Data -lGaspardData -lInteractionCoordinates \ - -lInitialConditions -lReaction -LDFLAGS += -L$(CLHEP_LIB_DIR) -l$(CLHEP_LIB) +LDFLAGS+= -L$(NPLIB)/lib -lVDetector -lIORoot -lReaction -lEnergyLoss \ + -lMust2Data -lMust2Physics \ + -lGaspardData -lGaspardPhysics \ + -lAnnularS1Data \ + -lInitialConditions -lInteractionCoordinates +LDFLAGS+= -L$(CLHEP_LIB_DIR) -l$(CLHEP_LIB) SRC= $(wildcard *.cc) -INC= $(wildcart $(NPAINCLUDES)/*.hh) OBJ=$(SRC:.cc=.o) #all:$(EXEC) # @$(CPP) -o $@ -c $< $(CXXFLAGS) -Analysis:$(OBJ) $(INC) +Analysis:$(OBJ) @$(CPP) $(CXXFLAGS) -o $@ $^ $(LDFLAGS) mv Analysis ../Analysis diff --git a/NPAnalysis/Gaspard/src/ObjectManager.hh b/NPAnalysis/Gaspard/src/ObjectManager.hh deleted file mode 100644 index 3e2bd0b1deb511b91eec39ccf0d7b8102d370d5d..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/src/ObjectManager.hh +++ /dev/null @@ -1,84 +0,0 @@ -// You can use this file to declare your spectra, file, energy loss , ... and whatever you want. -// This way you can remove all unnecessary declaration in the main programm. -// In order to help debugging and organising we use Name Space. - - -// -----------------------------------DOUBLE, INT, BOOL-------------------------------------------- - -namespace VARIABLE - { - double Xs,Ys,Zs,N ; - double Theta = 0 ; - double Ex = 0 ; - double ThetaCM = 0 ; - int X,Y ; - bool check_light = false ; - - } - -using namespace VARIABLE ; -// ---------------------------------------------------------------------------------------------- - - - -// -----------------------------------GRAPH------------------------------------------------------ -#include <TObject.h> -#include <TH1.h> -#include <TH1F.h> -#include <TH2.h> -#include <TH2F.h> -#include <TGraph2D.h> - -namespace GRAPH - { - // Declare your Spectra here: - const int NbElements = 40; - - TH2F *hXY1 = new TH2F("hXY1","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY2 = new TH2F("hXY2","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY3 = new TH2F("hXY3","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY4 = new TH2F("hXY4","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY5 = new TH2F("hXY5","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hXY6 = new TH2F("hXY6","", 128 , 1 , 128 , 128 , 1 , 128 ) ; - TH2F *hHitPattern[NbElements]; - - TH1F *hEnerX = new TH1F("hEnerX", "Energy X", 1000, -1, 10); - TH1F *hEnerY = new TH1F("hEnerY", "Energy Y", 1000, -1, 10); - TH2F *hEnerXY = new TH2F("hEnerXY", "Energy X vs Y", 100, -1, 10, 100, -1, 10); - TH2F *hEnerTheta = new TH2F("hEnerTheta", "Energy v.s. Theta", 90, 0, 180, 400, 0, 100); - - TH2F *hXY = new TH2F("hXY","", 300 , -130 , 130 , 300 , -130 , 130 ) ; - TH2F *hEA = new TH2F("hEA" ,"", 120 , 0 , 60 , 260 , 0 , 130 ) ; - TH1F *hEx = new TH1F("hEx", "E* (MeV)", 10000, -10, 40); - TH1F *hEHexa = new TH1F("hEHexa","EHexa (MeV)", 1000 , -5 , 5 ) ; - TH1F *hThetaCM = new TH1F("hThetaCM","ThetaCM (deg)", 30 , 0 , 60 ) ; - TH1F *hTheta = new TH1F("hTheta","Theta Lab (deg)",180,0,180) ; - TH2F *hThetaLabCM = new TH2F("hThetaLabCM","", 90 , 0 , 180 , 90 , 0 , 180 ) ; - - // Declare your Cut here: - - } - -using namespace GRAPH ; -// ---------------------------------------------------------------------------------------------- - - -// -----------------------------------ENERGY LOSS------------------------------------------------ -/* -#include "NPEnergyLossCorrection.h" -using namespace NPL ; -namespace ENERGYLOSS - { - -// EnergyLoss ProtonTarget = EnergyLoss ( "/home/Adrien/Desktop/dEdX.txt" , -// 100 , -// 1 ); - - } - -using namespace ENERGYLOSS ; -*/ -// ---------------------------------------------------------------------------------------------- - - - diff --git a/NPAnalysis/Gaspard/src/must.cc b/NPAnalysis/Gaspard/src/must.cc deleted file mode 100644 index 4ccdeb7ed9e1b749e14dde0e1b5e50ca7126ca28..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/src/must.cc +++ /dev/null @@ -1,245 +0,0 @@ -#include "must.hh" -#include <cmath> -#include <iostream> - -using namespace std; - -/////////////////////////////////////////////////////// - -//////// Constructor by Spherical Coordinate -must::must(double theta, double phi, double distance, double beta_u , double beta_v, double beta_w) - { - double Pi = 3.141592654 ; - - // convert from degree to radian: - theta = theta * Pi/180. ; - phi = phi * Pi/180. ; - - //Vector U on Telescope Face (paralelle to Y Strip) (NB: remember that Y strip are allong X axis) - Hep3Vector U ; - //Vector V on Telescope Face (parallele to X Strip) - Hep3Vector V ; - //Vector W normal to Telescope Face (pointing CsI) - Hep3Vector W ; - //Vector position of Telescope Face center - Hep3Vector C ; - - /*if(theta==180 && phi==90) - { - C = Hep3Vector (0,0,distance) ; - U = Hep3Vector (1,0,0) ; - V = Hep3Vector (0,1,0) ; - W = Hep3Vector (0,0,1) ; - }*/ - - if(theta==0 && phi==0) - { - C = Hep3Vector (0,0,distance) ; - U = Hep3Vector (1,0,0) ; - V = Hep3Vector (0,1,0) ; - W = Hep3Vector (0,0,1) ; - } - - else - { - C = Hep3Vector ( distance * sin(theta) * cos(phi) , - distance * sin(theta) * sin(phi) , - distance * cos(theta) ); - - Hep3Vector Y = Hep3Vector( cos(theta) * cos(phi), - cos(theta) * sin(phi), - -sin(theta)); - - W = C.unit() ; -// U = W.cross( Hep3Vector (0,1,0) ) ; - U = W.cross( Y ) ; - V = W.cross(U); - - U = U.unit(); - V = V.unit(); - - U.rotate( U , beta_u * Pi/180. ) ; - V.rotate( U , beta_u * Pi/180. ) ; - - U.rotate( V , beta_v * Pi/180. ) ; - V.rotate( V , beta_v * Pi/180. ) ; - - U.rotate( W , beta_w * Pi/180. ) ; - V.rotate( W , beta_w * Pi/180. ) ; - } - -// double Face = 98 ; //mm - double Face = 50 ; //mm - double NumberOfStrip = 128 ; - double StripPitch = Face/NumberOfStrip ; //mm - - vector<double> lineX ; vector<double> lineY ; vector<double> lineZ ; - double X , Y , Z ; - - //Moving C to the 1.1 corner: - C.setX( C.x() - ( Face/2 - StripPitch/2 ) * ( V.x() + U.x() ) ) ; - C.setY( C.y() - ( Face/2 - StripPitch/2 ) * ( V.y() + U.y() ) ) ; - C.setZ( C.z() - ( Face/2 - StripPitch/2 ) * ( V.z() + U.z() ) ) ; - - for( int i = 0 ; i < 128 ; i++ ) - { - - lineX.clear() ; - lineY.clear() ; - lineZ.clear() ; - - for( int j = 0 ; j < 128 ; j++ ) - { - - X = C.x() + StripPitch * ( U.x()*i + V.x()*j ) ; - Y = C.y() + StripPitch * ( U.y()*i + V.y()*j ) ; - Z = C.z() + StripPitch * ( U.z()*i + V.z()*j ) ; - - lineX.push_back(X) ; - lineY.push_back(Y) ; - lineZ.push_back(Z) ; - } - - StripPositionX.push_back(lineX) ; - StripPositionY.push_back(lineY) ; - StripPositionZ.push_back(lineZ) ; - - } - } - -/////////////////////////////////////////////////////// - -must::must( Hep3Vector C_X1_Y1 , - Hep3Vector C_X128_Y1 , - Hep3Vector C_X1_Y128 , - Hep3Vector C_X128_Y128 , - Hep3Vector TargetPosition = Hep3Vector(0,0,0) ) - { - - // Vector U on Telescope Face (paralelle to Y Strip) (NB: remember that Y strip are allong X axis) - Hep3Vector U = C_X128_Y1 - C_X1_Y1 ; - U = U.unit() ; - - // Vector V on Telescope Face (parallele to X Strip) - Hep3Vector V = C_X1_Y128 - C_X1_Y1 ; - V = V.unit() ; - - // Position Vector of Strip Center - Hep3Vector StripCenter ; - // Position Vector of X=1 Y=1 Strip - Hep3Vector Strip_1_1 ; - - // Geometry Parameter - double Face = 98 ; //mm - double NumberOfStrip = 128 ; - double StripPitch = Face/NumberOfStrip ; //mm - - // Buffer object to fill Position Array - vector<double> lineX ; vector<double> lineY ; vector<double> lineZ ; - - // Moving StripCenter to 1.1 corner: - Strip_1_1 = C_X1_Y1 + (U+V) * StripPitch/2 ; - - for( int i = 0 ; i < 128 ; i++ ) - { - lineX.clear() ; - lineY.clear() ; - lineZ.clear() ; - - for( int j = 0 ; j < 128 ; j++ ) - { - StripCenter = Strip_1_1 + StripPitch*( i*U + j*V ) ; - //StripCenter += -TargetPosition ; - - lineX.push_back( StripCenter.x() ) ; - lineY.push_back( StripCenter.y() ) ; - lineZ.push_back( StripCenter.z() ) ; - } - - StripPositionX.push_back(lineX) ; - StripPositionY.push_back(lineY) ; - StripPositionZ.push_back(lineZ) ; - - } - } - -/////////////////////////////////////////////////////// -double must::GetStripPositionX(int X , int Y) - { - return StripPositionX[X][Y]; - } -/////////////////////////////////////////////////////// -double must::GetStripPositionY(int X , int Y) - { - return StripPositionY[X][Y]; - } -/////////////////////////////////////////////////////// -double must::GetStripPositionZ(int X , int Y) - { - return StripPositionZ[X][Y]; - } -/////////////////////////////////////////////////////// -void must::Print() - { - - for( int i = 0 ; i < 128 ; i++ ) - { - - for( int j = 0 ; j < 128 ; j++ ) - { - cout << i+1 << " "<< j+1 << " " - << StripPositionX[i][j] << " " - << StripPositionY[i][j] << " " - << StripPositionZ[i][j] << " " - << endl ; - } - } - } - -/////////////////////////////////////////////////////// -///////////////// Class Array ///////////////////////// -/////////////////////////////////////////////////////// -array::array() -{NumberOfTelescope=0;} -/////////////////////////////////////////////////////// -void array::addTelescope( double theta , double phi , double distance , double beta_u , double beta_v, double beta_w) - { - must myMust(theta , phi , distance , beta_u , beta_v, beta_w) ; - NumberOfTelescope++ ; - myArray.push_back(myMust) ; - } - -/////////////// -void array::addTelescope(Hep3Vector A , Hep3Vector B , Hep3Vector C , Hep3Vector D ) - { - must myMust( A , B , C , D) ; - NumberOfTelescope++ ; - myArray.push_back(myMust) ; - } - - -/////////////////////////////////////////////////////// -double array::GetStripPositionX( int N , int X , int Y ) - { - return myArray[N-1].GetStripPositionX( X-1 , Y-1 ); - } - -/////////////////////////////////////////////////////// -double array::GetStripPositionY( int N , int X , int Y ) - { - return myArray[N-1].GetStripPositionY( X-1 , Y-1 ); - } - -/////////////////////////////////////////////////////// -double array::GetStripPositionZ( int N , int X , int Y ) - { - return myArray[N-1].GetStripPositionZ( X-1 , Y-1 ); - } - -/////////////////////////////////////////////////////// -void array::Print() - { - for ( int i = 0 ; i < myArray.size() ; i++ ) - { cout << " Telescope : " << i+1 << endl ; myArray.at(i).Print() ; } - } - diff --git a/NPAnalysis/Gaspard/src/must.hh b/NPAnalysis/Gaspard/src/must.hh deleted file mode 100644 index d3fa3a733fc029e7ed403f3a944e9dd68fc390c4..0000000000000000000000000000000000000000 --- a/NPAnalysis/Gaspard/src/must.hh +++ /dev/null @@ -1,75 +0,0 @@ -#ifndef must_h -#define must_h 1 - -#include <vector> -#include <TFile.h> -#include <TH2.h> -#include <TH2F.h> -#include "CLHEP/Vector/ThreeVector.h" -using namespace std ; -using namespace CLHEP ; - -class must -{ -public: - must( double theta , double phi , double distance , double beta_u , double beta_v, double beta_w) ; - must( Hep3Vector C_X1_Y1 , - Hep3Vector C_X128_Y1 , - Hep3Vector C_X1_Y128 , - Hep3Vector C_X128_Y128 , - Hep3Vector TargetPosition ); - - double GetStripPositionX( int X , int Y ) ; - double GetStripPositionY( int X , int Y ) ; - double GetStripPositionZ( int X , int Y ) ; - void Print() ; - -private: - - vector< vector < double > > StripPositionX ; - vector< vector < double > > StripPositionY ; - vector< vector < double > > StripPositionZ ; - -}; - - -class array -{ -public: - array(); - void addTelescope( double theta , double phi , double distance , double beta_u , double beta_v , double beta_w) ; - void addTelescope( Hep3Vector C_X1_Y1 , - Hep3Vector C_X128_Y1 , - Hep3Vector C_X1_Y128 , - Hep3Vector C_X128_Y128 ); - - double GetStripPositionX( int N , int X , int Y ) ; - double GetStripPositionY( int N , int X , int Y ) ; - double GetStripPositionZ( int N , int X , int Y ) ; - double GetNumberOfTelescope() {return NumberOfTelescope;} ; - void Print() ; - -private: - vector< must > myArray ; - int NumberOfTelescope ; - -public: // Histogram - -/*// Hit Density - vector<TH2F*> HitDensity ; - - // Time of Flight Identification - vector<TH2F*> TOF ; - - // EDE Identification - vector<TH2F*> EDE ; - - // Kinematic Line - KinematicLine ; - - // Excitation Energy - ExcitationEnergy ;*/ - -}; - -#endif diff --git a/NPAnalysis/Gaspard/macros/ControlSimu.C b/NPAnalysis/macros/ControlSimu.C similarity index 70% rename from NPAnalysis/Gaspard/macros/ControlSimu.C rename to NPAnalysis/macros/ControlSimu.C index 56f2ee46029f7451fcf69349a615f5cd3b5a17bc..bb7b8ddcbe521156b05f6c42264ec1922d0e8456 100644 --- a/NPAnalysis/Gaspard/macros/ControlSimu.C +++ b/NPAnalysis/macros/ControlSimu.C @@ -1,3 +1,30 @@ +/***************************************************************************** + * Copyright (C) 2009 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: N. de Sereville contact address: deserevi@ipno.in2p3.fr * + * * + * Creation Date : 22/07/09 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * + This macro displays everything concerning the incident beam and the * + * emitted particle from NPSimulation * + * * + * + Use in a ROOT session: * + * .x ControlSimu.C("FileToAnalyse") * + * * + * * + *---------------------------------------------------------------------------* + * Comment: * + * * + * * + *****************************************************************************/ + #include <iostream> #include "TROOT.h" diff --git a/NPAnalysis/Gaspard/macros/GeometricalEfficiency.C b/NPAnalysis/macros/GeometricalEfficiency.C similarity index 50% rename from NPAnalysis/Gaspard/macros/GeometricalEfficiency.C rename to NPAnalysis/macros/GeometricalEfficiency.C index 90ca0e83b863b7b45116682d97a9fafb246068af..5bbc6f5864000c292f58a2741c5f69160e0ed90f 100644 --- a/NPAnalysis/Gaspard/macros/GeometricalEfficiency.C +++ b/NPAnalysis/macros/GeometricalEfficiency.C @@ -1,3 +1,31 @@ +/***************************************************************************** + * Copyright (C) 2009 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: N. de Sereville contact address: deserevi@ipno.in2p3.fr * + * * + * Creation Date : 22/07/09 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * + This macro calculates the geometrical efficiency for a given geometry* + * NPSimulation should have been run with an isotropic source event * + * generator. * + * * + * + Use in a ROOT session: * + * .x GeometricalEfficiency.C("FileToAnalyse") * + * * + * * + *---------------------------------------------------------------------------* + * Comment: * + * * + * * + *****************************************************************************/ + #include <iostream> #include "TROOT.h" diff --git a/NPAnalysis/macros/README b/NPAnalysis/macros/README new file mode 100644 index 0000000000000000000000000000000000000000..de21ad1250c784e67d5ae2eadfa26e9c4ec9841f --- /dev/null +++ b/NPAnalysis/macros/README @@ -0,0 +1,24 @@ +/***************************************************************************** + * Copyright (C) 2009 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: N. de Sereville contact address: deserevi@ipno.in2p3.fr * + * * + * Creation Date : 11/09/09 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * + The ROOT macros found in this directory are independent of the * + * detector set-up and can be run for any detector. * + * * + * * + *---------------------------------------------------------------------------* + * Comment: * + * * + * * + *****************************************************************************/ + diff --git a/NPAnalysis/Gaspard/macros/TestAngleReconstruction.C b/NPAnalysis/macros/TestAngleReconstruction.C similarity index 68% rename from NPAnalysis/Gaspard/macros/TestAngleReconstruction.C rename to NPAnalysis/macros/TestAngleReconstruction.C index 217efe959dae85053ab602d6efdc92ff9b05df2b..b29fc482f77dca795824ee3e2e2ad710f59d939d 100644 --- a/NPAnalysis/Gaspard/macros/TestAngleReconstruction.C +++ b/NPAnalysis/macros/TestAngleReconstruction.C @@ -1,3 +1,31 @@ +/***************************************************************************** + * Copyright (C) 2009 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: N. de Sereville contact address: deserevi@ipno.in2p3.fr * + * * + * Creation Date : 22/07/09 * + * Last update : * + *---------------------------------------------------------------------------* + * Decription: * + * + This macro studues the impact of the beam emittance on the angle * + * determination. * + * * + * + Use in a ROOT session: * + * .x TestAngleReconstruction.C("FileToAnalyse") * + * * + * * + *---------------------------------------------------------------------------* + * Comment: * + * + This macro has not been fully tested yet. Be careful * + * * + * * + *****************************************************************************/ + #include <iostream> #include "TROOT.h" @@ -11,7 +39,6 @@ #include "TCanvas.h" #include "TMath.h" -#include "TGaspardTrackerData.h" #include "TInitialConditions.h" #include "TInteractionCoordinates.h" @@ -23,13 +50,9 @@ void TestAngleReconstruction(const char * fname = "mySimul") TString inFileName = fname; if (!inFileName.Contains("root")) inFileName += ".root"; TFile *inFile = new TFile(path + inFileName); - TTree *tree = (TTree*) inFile->Get("EventTree"); + TTree *tree = (TTree*) inFile->Get("SimulatedTree"); // Connect the branches of the TTree and activate then if necessary - // TGaspardTrackerData branch - TGaspardTrackerData *gpdTrkData = 0; - tree->SetBranchAddress("GASPARD", &gpdTrkData); - tree->SetBranchStatus("GASPARD", 0); // TInitialConditions branch TInitialConditions *initCond = 0; tree->SetBranchAddress("InitialConditions", &initCond); diff --git a/NPLib/GASPARD/GaspardTracker.cxx b/NPLib/GASPARD/GaspardTracker.cxx index d9971e51b4cbcb4c8a5b24d90a3be91dc2b0131a..c6afac495de05ac57b9fd2a2dc6e1f795bac21c0 100644 --- a/NPLib/GASPARD/GaspardTracker.cxx +++ b/NPLib/GASPARD/GaspardTracker.cxx @@ -244,9 +244,9 @@ void GaspardTracker::ReadConfiguration(string Path) // with angle method else if ( check_Theta && check_Phi && check_R && check_beta ) { - AddModuleSquare(R , - Theta , + AddModuleSquare(Theta , Phi , + R , beta_u , beta_v , beta_w ); @@ -398,12 +398,12 @@ void GaspardTracker::ReadConfiguration(string Path) // with angle method else if ( check_Theta && check_Phi && check_R && check_beta ) { - AddModuleDummyShape(R , - Theta , - Phi , - beta_u , - beta_v , - beta_w ); + AddModuleDummyShape(Theta, + Phi, + R, + beta_u, + beta_v, + beta_w); } // reset boolean flag for point positioning diff --git a/NPLib/GASPARD/GaspardTracker.h b/NPLib/GASPARD/GaspardTracker.h index 7ffc4f8db5a6098750ffcfe13f7b7d63f9821382..e0d25a08bab1e7316a715f67c2a9dba4fd9df4d3 100644 --- a/NPLib/GASPARD/GaspardTracker.h +++ b/NPLib/GASPARD/GaspardTracker.h @@ -106,14 +106,14 @@ public: double beta_w); // Getters to retrieve the (X,Y,Z) coordinates of a pixel defined by strips (X,Y) - double GetStripPositionX(int N ,int X ,int Y) { return m_StripPositionX[N-1][X-1][Y-1]; }; - double GetStripPositionY(int N ,int X ,int Y) { return m_StripPositionY[N-1][X-1][Y-1]; }; - double GetStripPositionZ(int N ,int X ,int Y) { return m_StripPositionZ[N-1][X-1][Y-1]; }; - double GetNumberOfModule() { return m_NumberOfModule; }; + double GetStripPositionX(int N ,int X ,int Y) { return m_StripPositionX[N-1][X-1][Y-1]; } + double GetStripPositionY(int N ,int X ,int Y) { return m_StripPositionY[N-1][X-1][Y-1]; } + double GetStripPositionZ(int N ,int X ,int Y) { return m_StripPositionZ[N-1][X-1][Y-1]; } + double GetNumberOfModule() { return m_NumberOfModule; } // Get Root input and output objects - TGaspardTrackerData* GetEventData() {return m_EventData;}; - TGaspardTrackerPhysics* GetEventPhysics() {return m_EventPhysics;}; + TGaspardTrackerData* GetEventData() {return m_EventData;} + TGaspardTrackerPhysics* GetEventPhysics() {return m_EventPhysics;} // To be called after a build Physical Event double GetEnergyDeposit(); diff --git a/NPLib/GASPARD/TGaspardTrackerPhysics.cxx b/NPLib/GASPARD/TGaspardTrackerPhysics.cxx index 6bca92aea82e27c4d075f453a5f9370a91a82093..51f79ca0e9fd8402969bda42adedccb79d4e6e10 100644 --- a/NPLib/GASPARD/TGaspardTrackerPhysics.cxx +++ b/NPLib/GASPARD/TGaspardTrackerPhysics.cxx @@ -136,7 +136,7 @@ void TGaspardTrackerPhysics::BuildPhysicalEvent(TGaspardTrackerData* Data) cout << "Warning: multiplicity in third stage greater than in firststage" << endl; } - // Analysis code here. + // Fill total energy TotalEnergy.push_back(EnergyTot); } else { diff --git a/NPLib/InteractionCoordinates/TInteractionCoordinates.cxx b/NPLib/InteractionCoordinates/TInteractionCoordinates.cxx index ac06880b847bbc456cc51c60a6d9a7758e33894f..25e71b4f36ec59f741a4089304621bf14fcd12b9 100644 --- a/NPLib/InteractionCoordinates/TInteractionCoordinates.cxx +++ b/NPLib/InteractionCoordinates/TInteractionCoordinates.cxx @@ -54,9 +54,9 @@ void TInteractionCoordinates::Clear() void TInteractionCoordinates::Dump() { - cout << "XXXXXXXXXXXXX Initial conditions XXXXXXXXXXXXXXXX" << endl; + cout << "XXXXXXXXXXXXX Interaction coordinates XXXXXXXXXXXXXXXX" << endl; - cout << "Vertex position : " << endl; + cout << "Interaction position : " << endl; cout << "\tX : " << fDetected_Position_X[0] << endl; cout << "\tY : " << fDetected_Position_Y[0] << endl; cout << "\tZ : " << fDetected_Position_Z[0] << endl;