From aae167653e81726e41e1e26bd5353aea96db4433 Mon Sep 17 00:00:00 2001
From: Pierre Morfouace <pierre.morfouace2@cea.fr>
Date: Tue, 26 Oct 2021 15:41:18 +0200
Subject: [PATCH] * Updating s455 project
---
Projects/s455/Analysis.cxx | 503 ++++++++++--------
Projects/s455/Analysis.h | 10 +
.../calibration/SofSci/SofSci_physics.cal | 4 +-
.../SofTrim/SofTrim_SectionAlign.cal | 4 +-
Projects/s455/configs/ConfigSofSci.dat | 4 +-
Projects/s455/macro/GetChargeDist.C | 352 ++++++++++++
6 files changed, 648 insertions(+), 229 deletions(-)
create mode 100644 Projects/s455/macro/GetChargeDist.C
diff --git a/Projects/s455/Analysis.cxx b/Projects/s455/Analysis.cxx
index 9b0c02b6b..f60e2823c 100644
--- a/Projects/s455/Analysis.cxx
+++ b/Projects/s455/Analysis.cxx
@@ -56,6 +56,7 @@ void Analysis::Init(){
void Analysis::TreatEvent(){
ReInitValue();
//cout << "************" << endl;
+ RunID = fRunID;
BeamAnalysis();
unsigned int sofsci_size = SofSci->DetectorNbr.size();
@@ -67,8 +68,6 @@ void Analysis::TreatEvent(){
FissionFragmentAnalysis();
}
-
-
}
////////////////////////////////////////////////////////////////////////////////
@@ -127,7 +126,6 @@ void Analysis::FissionFragmentAnalysis(){
int multR = mult3 + mult4;
if(softwim_size>1){
- //if( (multL>15 && multR>15) || (multR==32) || (multL==32)){
if( (mult1>1 && mult1<17) || (mult2>1 && mult2<17) || (mult3>1 && mult3<17) || (mult4>1 && mult4<17)){
for(unsigned int i=0; i< softwim_size; i++){
int sec = SofTwim->SectionNbr[i];
@@ -247,9 +245,6 @@ void Analysis::FissionFragmentAnalysis(){
// Z calibration //
- double p0 = 2.80063;
- double p1 = 6.91985e-2;
- double p2 = 1.01598e-7;
double Z1=-1;
double Z2=-1;
double Zsum=-1;
@@ -281,9 +276,9 @@ void Analysis::FissionFragmentAnalysis(){
}
if(Z1>0 && Z2>0){
- Z1 = p0 + p1*Z1 + p2*Z1*Z1;
- Z2 = p0 + p1*Z2 + p2*Z2*Z2;
-
+ Z1 = fZff_p0 + fZff_p1*Z1 + fZff_p2*Z1*Z1;
+ Z2 = fZff_p0 + fZff_p1*Z2 + fZff_p2*Z2*Z2;
+
Z1 = sqrt(Z1);
Z2 = sqrt(Z2);
@@ -301,15 +296,6 @@ void Analysis::FissionFragmentAnalysis(){
SofFF->SetBeta(Beta3);
SofFF->SetBeta(Beta4);
- /*SofFF->SetZ(E1);
- SofFF->SetZ(E2);
- SofFF->SetZ(E3);
- SofFF->SetZ(E4);*/
-
- /*SofFF->SetZ( sqrt(p0+p1*E1+p2*E1*E1) );
- SofFF->SetZ( sqrt(p0+p1*E2+p2*E2*E2) );
- SofFF->SetZ( sqrt(p0+p1*E3+p2*E3*E3) );
- SofFF->SetZ( sqrt(p0+p1*E4+p2*E4*E4) );*/
SofFF->SetZ(Z1);
SofFF->SetZ(Z2);
@@ -322,234 +308,305 @@ void Analysis::FissionFragmentAnalysis(){
SofFF->SetIntZsum(iZsum);
}
}
- }
+}
- ////////////////////////////////////////////////////////////////////////////////
- void Analysis::BeamAnalysis(){
- unsigned int sofsci_size = SofSci->DetectorNbr.size();
- double Y_p0 = 23943.8;
- double Y_p1 = 12.362;
-
- double alpha = 2.09862;
- double beta = 0.972689;
- double Z_p0 = alpha+beta*-9.31873;
- double Z_p1 = beta*0.564459;
-
- if(sofsci_size==2){
- double beta = SofSci->Beta[0];
- //cout << "Set beta to " << beta << endl;
- SofTrim->SetBeta(beta);
- SofTrim->BuildSimplePhysicalEvent();
- double Zbeam,Qmax,Theta;
- if(SofTrim->EnergySection.size()>0){
- Zbeam = SofTrim->GetMaxEnergySection();
- Qmax = DetermineQmax();
- Theta = SofTrim->Theta[0];
-
- double TofFromS2 = SofSci->CalTof[0];
- double velocity_mns = SofSci->VelocityMNs[0];
- double Beta = SofSci->Beta[0];
- double XS2 = SofSci->PosMm[0];
- //double XCC = SofSci->PosMm[1];
- double XCC=0;
- double YCC=0;
- for(unsigned int i=0; i<SofMwpc->DetectorNbr.size(); i++){
- if(SofMwpc->DetectorNbr[i]==1){
- XCC = SofMwpc->PositionX1[i];
- YCC = SofMwpc->PositionY[i];
- }
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::BeamAnalysis(){
+ unsigned int sofsci_size = SofSci->DetectorNbr.size();
+ if(sofsci_size==2){
+ double beta = SofSci->Beta[0];
+ //cout << "Set beta to " << beta << endl;
+ SofTrim->SetBeta(beta);
+ SofTrim->BuildSimplePhysicalEvent();
+ double Zbeam,Qmax,Theta;
+ if(SofTrim->EnergySection.size()>0){
+ Zbeam = SofTrim->GetMaxEnergySection();
+ Qmax = DetermineQmax();
+ Theta = SofTrim->Theta[0];
+
+ double TofFromS2 = SofSci->CalTof[0];
+ double velocity_mns = SofSci->VelocityMNs[0];
+ double Beta = SofSci->Beta[0];
+ double XS2 = SofSci->PosMm[0];
+ //double XCC = SofSci->PosMm[1];
+ double XCC=0;
+ double YCC=0;
+ for(unsigned int i=0; i<SofMwpc->DetectorNbr.size(); i++){
+ if(SofMwpc->DetectorNbr[i]==1){
+ XCC = SofMwpc->PositionX1[i];
+ YCC = SofMwpc->PositionY[i];
}
-
- double LS2;
- LS2 = fLS2_0;//*(1 + fK_LS2*Theta);
- velocity_mns = LS2/TofFromS2;
- Beta = velocity_mns * m/ns / NPUNITS::c_light;
- double Gamma = 1./(TMath::Sqrt(1 - TMath::Power(Beta,2)));
- double Brho = fBrho0 * (1 - XS2/fDS2 - XCC/fDCC);
- double AoQ = Brho / (3.10716*Gamma*Beta);
- double A = AoQ * Qmax;
-
- // Y dependence correction //
- Zbeam = Zbeam/(Y_p0 + Y_p1*YCC)*Y_p0;
-
- // Z calibration //
- Zbeam = Z_p0 + Z_p1*sqrt(Zbeam);
- //Zbeam = 2.09862 + 0.972689*Zbeam;
-
- // Filling Beam tree
- SofBeamID->SetZbeam(Zbeam);
- SofBeamID->SetQmax(rand.Gaus(Qmax,0.15));
- SofBeamID->SetAoQ(AoQ);
- SofBeamID->SetAbeam(A);
- SofBeamID->SetBeta(Beta);
- SofBeamID->SetGamma(Gamma);
- SofBeamID->SetBrho(Brho);
- SofBeamID->SetXS2(XS2);
- SofBeamID->SetXCC(XCC);
- SofBeamID->SetYCC(YCC);
}
- }
- }
- ////////////////////////////////////////////////////////////////////////////////
- void Analysis::LoadCut(){
- TString input_path = "./calibration/SofTrim/cut/";
-
- TString rootfile;
- TString cutfile;
- TFile* file;
- for(int i=0; i<3; i++){
- // Q=77
- rootfile = Form("cutsec%iQ77.root",i+1);
- cutfile = input_path + rootfile;
- file = new TFile(cutfile);
- cutQ77[i] = (TCutG*) file->Get(Form("cutsec%iQ77",i+1));
- // Q=78
- rootfile = Form("cutsec%iQ78.root",i+1);
- cutfile = input_path + rootfile;
- file = new TFile(cutfile);
- cutQ78[i] = (TCutG*) file->Get(Form("cutsec%iQ78",i+1));
- // Q=79
- rootfile = Form("cutsec%iQ79.root",i+1);
- cutfile = input_path + rootfile;
- file = new TFile(cutfile);
- cutQ79[i] = (TCutG*) file->Get(Form("cutsec%iQ79",i+1));
- // Q=80
- rootfile = Form("cutsec%iQ80.root",i+1);
- cutfile = input_path + rootfile;
- file = new TFile(cutfile);
- cutQ80[i] = (TCutG*) file->Get(Form("cutsec%iQ80",i+1));
+ double LS2;
+ LS2 = fLS2_0;//*(1 + fK_LS2*Theta);
+ velocity_mns = LS2/TofFromS2;
+ Beta = velocity_mns * m/ns / NPUNITS::c_light;
+ double Gamma = 1./(TMath::Sqrt(1 - TMath::Power(Beta,2)));
+ double Brho = fBrho0 * (1 - XS2/fDS2 - XCC/fDCC);
+ double AoQ = Brho / (3.10716*Gamma*Beta);
+ double A = AoQ * Qmax;
+
+ // Y dependence correction //
+ double Y_p0 = 23943.8;
+ double Y_p1 = 12.362;
+ Zbeam = Zbeam/(Y_p0 + Y_p1*YCC)*Y_p0;
+
+ // Z calibration //
+ Zbeam = fZbeam_p0 + fZbeam_p1*Zbeam + fZbeam_p2*Zbeam*Zbeam;
+ Zbeam = sqrt(Zbeam);
+
+ // Last beta correction //
+ double Beta_norm = 0.8355;
+ Zbeam = Zbeam/(fZBeta_p0 + fZBeta_p1*Beta)*(fZBeta_p0 + fZBeta_p1*Beta_norm);
+
+ // Filling Beam tree
+ SofBeamID->SetZbeam(Zbeam);
+ SofBeamID->SetQmax(rand.Gaus(Qmax,0.15));
+ SofBeamID->SetAoQ(AoQ);
+ SofBeamID->SetAbeam(A);
+ SofBeamID->SetBeta(Beta);
+ SofBeamID->SetGamma(Gamma);
+ SofBeamID->SetBrho(Brho);
+ SofBeamID->SetXS2(XS2);
+ SofBeamID->SetXCC(XCC);
+ SofBeamID->SetYCC(YCC);
}
}
+}
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::LoadCut(){
+ TString input_path = "./calibration/SofTrim/cut/";
+
+ TString rootfile;
+ TString cutfile;
+ TFile* file;
+ for(int i=0; i<3; i++){
+ // Q=77
+ rootfile = Form("cutsec%iQ77.root",i+1);
+ cutfile = input_path + rootfile;
+ file = new TFile(cutfile);
+ cutQ77[i] = (TCutG*) file->Get(Form("cutsec%iQ77",i+1));
+ // Q=78
+ rootfile = Form("cutsec%iQ78.root",i+1);
+ cutfile = input_path + rootfile;
+ file = new TFile(cutfile);
+ cutQ78[i] = (TCutG*) file->Get(Form("cutsec%iQ78",i+1));
+ // Q=79
+ rootfile = Form("cutsec%iQ79.root",i+1);
+ cutfile = input_path + rootfile;
+ file = new TFile(cutfile);
+ cutQ79[i] = (TCutG*) file->Get(Form("cutsec%iQ79",i+1));
+ // Q=80
+ rootfile = Form("cutsec%iQ80.root",i+1);
+ cutfile = input_path + rootfile;
+ file = new TFile(cutfile);
+ cutQ80[i] = (TCutG*) file->Get(Form("cutsec%iQ80",i+1));
+ }
+}
- ////////////////////////////////////////////////////////////////////////////////
- void Analysis::LoadSpline(){
- TString input_path = "./calibration/SofTwim/spline/";
-
- TString rootfile = input_path + "spline_beta.root";
- TFile* ifile = new TFile(rootfile,"read");
- TString splinename;
- if(ifile->IsOpen()){
- cout << "Loading Beta spline for fission fragment analysis..." << endl;
- for(int i=0; i<4; i++){
- splinename = Form("spline_beta_sec%i",i+1);
- fcorr_z_beta[i] = (TSpline3*) ifile->FindObjectAny(splinename);
- }
- ifile->Close();
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::LoadSpline(){
+ TString input_path = "./calibration/SofTwim/spline/";
+
+ TString rootfile = input_path + "spline_beta.root";
+ TFile* ifile = new TFile(rootfile,"read");
+
+ TString splinename;
+ if(ifile->IsOpen()){
+ cout << "Loading Beta spline for fission fragment analysis..." << endl;
+ for(int i=0; i<4; i++){
+ splinename = Form("spline_beta_sec%i",i+1);
+ fcorr_z_beta[i] = (TSpline3*) ifile->FindObjectAny(splinename);
}
- else
- cout << "File " << rootfile << " not found!" << endl;
-
- //*** ***//
- rootfile = input_path + "spline_dt.root";
- ifile = new TFile(rootfile,"read");
-
- if(ifile->IsOpen()){
- cout << "Loading DT spline for fission fragment analysis..." << endl;
- for(int i=0; i<4; i++){
- splinename = Form("spline_dt_sec%i",i+1);
- fcorr_z_dt[i] = (TSpline3*) ifile->FindObjectAny(splinename);
- }
- ifile->Close();
+ ifile->Close();
+ }
+ else
+ cout << "File " << rootfile << " not found!" << endl;
+
+ //*** ***//
+ rootfile = input_path + "spline_dt.root";
+ ifile = new TFile(rootfile,"read");
+
+ if(ifile->IsOpen()){
+ cout << "Loading DT spline for fission fragment analysis..." << endl;
+ for(int i=0; i<4; i++){
+ splinename = Form("spline_dt_sec%i",i+1);
+ fcorr_z_dt[i] = (TSpline3*) ifile->FindObjectAny(splinename);
}
- else
- cout << "File " << rootfile << " not found!" << endl;
+ ifile->Close();
+ }
+ else
+ cout << "File " << rootfile << " not found!" << endl;
+
+}
+////////////////////////////////////////////////////////////////////////////////
+int Analysis::DetermineQmax(){
+ int Qmax;
+ int Qsec[3];
+
+ unsigned int size = SofTrim->EnergySection.size();
+ for(unsigned int i=0; i<size; i++){
+ int SecNbr = SofTrim->SectionNbr[i];
+ double Theta = SofTrim->Theta[i];
+ double Esec = SofTrim->EnergySection[i];
+
+
+ if(cutQ77[SecNbr-1]->IsInside(Theta,Esec))
+ Qsec[SecNbr-1] = 77;
+ else if(cutQ78[SecNbr-1]->IsInside(Theta,Esec))
+ Qsec[SecNbr-1] = 78;
+ else if(cutQ79[SecNbr-1]->IsInside(Theta,Esec))
+ Qsec[SecNbr-1] = 79;
+ else if(cutQ80[SecNbr-1]->IsInside(Theta,Esec))
+ Qsec[SecNbr-1] = 80;
+ //else if(cutQ81[SecNbr-1]->IsInside(Theta,Esec))
+ //Qsec[SecNbr-1] = 81;
}
- ////////////////////////////////////////////////////////////////////////////////
- int Analysis::DetermineQmax(){
- int Qmax;
- int Qsec[3];
-
- unsigned int size = SofTrim->EnergySection.size();
- for(unsigned int i=0; i<size; i++){
- int SecNbr = SofTrim->SectionNbr[i];
- double Theta = SofTrim->Theta[i];
- double Esec = SofTrim->EnergySection[i];
-
-
- if(cutQ77[SecNbr-1]->IsInside(Theta,Esec))
- Qsec[SecNbr-1] = 77;
- else if(cutQ78[SecNbr-1]->IsInside(Theta,Esec))
- Qsec[SecNbr-1] = 78;
- else if(cutQ79[SecNbr-1]->IsInside(Theta,Esec))
- Qsec[SecNbr-1] = 79;
- else if(cutQ80[SecNbr-1]->IsInside(Theta,Esec))
- Qsec[SecNbr-1] = 80;
- //else if(cutQ81[SecNbr-1]->IsInside(Theta,Esec))
- //Qsec[SecNbr-1] = 81;
- }
+ Qmax = max(Qsec[0],Qsec[1]);
+ Qmax = max(Qsec[2],Qmax);
- Qmax = max(Qsec[0],Qsec[1]);
- Qmax = max(Qsec[2],Qmax);
+ return Qmax;
+}
- return Qmax;
- }
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::End(){
+}
- ////////////////////////////////////////////////////////////////////////////////
- void Analysis::End(){
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::InitParameter(){
+ fLS2_0 = 135.614;
+ fDS2 = 8000;
+ fDCC = -10000;
+ fK_LS2 = -30e-8;
+
+ fRunID = 8;
+
+ // Beam parameter //
+ fZBeta_p0 = 1;
+ fZBeta_p1 = 0;
+
+ // FF paramter //
+ fZff_p0 = 2.80063;
+ fZff_p1 = 6.91985e-2;
+ fZff_p2 = 1.01598e-7;
+
+ if(fRunID==1 || fRunID==2){
+ //fBrho0 = 10.6813; // 180Hg
+ fBrho0 = 10.6955; // 180Hg
+ fZbeam_p0 = -5303.06;
+ fZbeam_p1 = 0.674945;
+ fZbeam_p2 = -8.32085e-6;
+
+ fZBeta_p0 = 2.73;
+ fZBeta_p1 = 74.7153;
}
-
- ////////////////////////////////////////////////////////////////////////////////
- void Analysis::InitParameter(){
- fLS2_0 = 135.614;
- fDS2 = 8000;
- fDCC = -10000;
- fK_LS2 = -30e-8;
-
-
- //fBrho0 = 14.1008; // 238U run 367
- //fBrho0 = 12.8719; // 238U run 368
- //fBrho0 = 12.3255; // 238U run 369
- //fBrho0 = 12.3352; // 238U run 422
- //fBrho0 = 10.9476; // 189Pb
- //fBrho0 = 10.9558; // 184Hg
- //fBrho0 = 10.8183; // 182Hg
- //fBrho0 = 10.6814; // 180Hg
- //fBrho0 = 10.8138; // 187Pb
- //fBrho0 = 11.3418; // 216Th
+ if(fRunID==3){
+ fBrho0 = 10.8183; // 182Hg
+ fZbeam_p0 = -2737.25;
+ fZbeam_p1 = 0.452017;
+ fZbeam_p2 = -3.48831e-6;
+ }
+ if(fRunID==4){
+ fBrho0 = 10.9558; // 184Hg
+ fZbeam_p0 = -5044.61;
+ fZbeam_p1 = 0.639986;
+ fZbeam_p2 = -7.3077e-6;
+ }
+ if(fRunID==5){
+ fBrho0 = 10.8138; // 187Pb
+ fZbeam_p0 = -2858.72;
+ fZbeam_p1 = 0.454064;
+ fZbeam_p2 = -3.36443e-6;
+ }
+ if(fRunID==6){
+ fBrho0 = 10.9476; // 189Pb
+ fZbeam_p0 = 1590.66;
+ fZbeam_p1 = 0.0956455;
+ fZbeam_p2 = 3.84585e-6;
+ }
+ if(fRunID==7){
+ fBrho0 = 10.6814; // 175Pt
+ fZbeam_p0 = 459.68;
+ fZbeam_p1 = 0.162277;
+ fZbeam_p2 = 3.10164e-6;
+ }
+ if(fRunID==8){
fBrho0 = 11.0864; // 204Fr
- //fBrho0 = 11.2712; // 207Fr
- //fBrho0 = 10.6814; // 175Pt
- //fBrho0 = 11.5067; // 221Pa
- //fBrho0 = 11.0955; // 199At run 423 & 424
- //fBrho0 = 10.9970; // 199At run 425 & 426
- //fBrho0 = 10.8697; //197At
+ fZbeam_p0 = 4122.94;
+ fZbeam_p1 = -0.119867;
+ fZbeam_p2 = 8.29115e-6;
+
+ fZBeta_p0 = 63.9575;
+ fZBeta_p1 = 25.1988;
}
-
- ////////////////////////////////////////////////////////////////////////////////
- void Analysis::ReInitValue(){
- SofBeamID->Clear();
- SofFF->Clear();
+ if(fRunID==9){
+ fBrho0 = 11.2712; // 207Fr
+ fZbeam_p0 = -1752.27;
+ fZbeam_p1 = 0.346018;
+ fZbeam_p2 = -8.64673e-7;
}
- ////////////////////////////////////////////////////////////////////////////////
- void Analysis::InitOutputBranch(){
- //RootOutput::getInstance()->GetTree()->Branch("Zbeam",&Zbeam,"Zbeam/D");
- RootOutput::getInstance()->GetTree()->Branch("SofBeamID","TSofBeamID",&SofBeamID);
- RootOutput::getInstance()->GetTree()->Branch("SofFissionFragment","TSofFissionFragment",&SofFF);
-
+ if(fRunID==10){
+ fBrho0 = 11.0955; // 199At run 423 & 424
+ fZbeam_p0 = -116.425;
+ fZbeam_p1 = 0.218256;
+ fZbeam_p2 = 1.62399e-6;
+ }
+ if(fRunID==11){
+ fBrho0 = 10.9970; // 199At run 425 & 426
+ fZbeam_p0 = -116.425;
+ fZbeam_p1 = 0.218256;
+ fZbeam_p2 = 1.62399e-6;
}
- ////////////////////////////////////////////////////////////////////////////////
- // Construct Method to be pass to the DetectorFactory //
- ////////////////////////////////////////////////////////////////////////////////
- NPL::VAnalysis* Analysis::Construct(){
- return (NPL::VAnalysis*) new Analysis();
+ if(fRunID==12){
+ fBrho0 = 10.8697; //197At
+ fZbeam_p0 = -2683.52;
+ fZbeam_p1 = 0.422551;
+ fZbeam_p2 = -2.44471e-6;
}
+ if(fRunID==13){
+ fBrho0 = 11.3418; // 216Th
+ }
+ if(fRunID==14){
+ fBrho0 = 11.5067; // 221Pa
+ }
+}
- ////////////////////////////////////////////////////////////////////////////////
- // Registering the construct method to the factory //
- ////////////////////////////////////////////////////////////////////////////////
- extern "C"{
- class proxy{
- public:
- proxy(){
- NPL::AnalysisFactory::getInstance()->SetConstructor(Analysis::Construct);
- }
- };
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::ReInitValue(){
+ SofBeamID->Clear();
+ SofFF->Clear();
+}
+////////////////////////////////////////////////////////////////////////////////
+void Analysis::InitOutputBranch(){
+ RootOutput::getInstance()->GetTree()->Branch("RunID",&RunID,"RunID/I");
+ RootOutput::getInstance()->GetTree()->Branch("SofBeamID","TSofBeamID",&SofBeamID);
+ RootOutput::getInstance()->GetTree()->Branch("SofFissionFragment","TSofFissionFragment",&SofFF);
+
+}
+////////////////////////////////////////////////////////////////////////////////
+// Construct Method to be pass to the DetectorFactory //
+////////////////////////////////////////////////////////////////////////////////
+NPL::VAnalysis* Analysis::Construct(){
+ return (NPL::VAnalysis*) new Analysis();
+}
+
+////////////////////////////////////////////////////////////////////////////////
+// Registering the construct method to the factory //
+////////////////////////////////////////////////////////////////////////////////
+extern "C"{
+ class proxy{
+ public:
+ proxy(){
+ NPL::AnalysisFactory::getInstance()->SetConstructor(Analysis::Construct);
+ }
+ };
+
+ proxy p;
+}
- proxy p;
- }
diff --git a/Projects/s455/Analysis.h b/Projects/s455/Analysis.h
index 447b56628..d4f2707a7 100644
--- a/Projects/s455/Analysis.h
+++ b/Projects/s455/Analysis.h
@@ -67,13 +67,23 @@ class Analysis: public NPL::VAnalysis{
TSofTofWPhysics* SofTofW;
TSofBeamID* SofBeamID;
TSofFissionFragment* SofFF;
+ int RunID;
private:
+ int fRunID;
double fLS2_0;
double fBrho0;
double fDS2;
double fDCC;
double fK_LS2;
+ double fZbeam_p0;
+ double fZbeam_p1;
+ double fZbeam_p2;
+ double fZff_p0;
+ double fZff_p1;
+ double fZff_p2;
+ double fZBeta_p0;
+ double fZBeta_p1;
TRandom3 rand;
TCutG* cutQ77[3];
diff --git a/Projects/s455/calibration/SofSci/SofSci_physics.cal b/Projects/s455/calibration/SofSci/SofSci_physics.cal
index b5d573448..be44a7ac2 100644
--- a/Projects/s455/calibration/SofSci/SofSci_physics.cal
+++ b/Projects/s455/calibration/SofSci/SofSci_physics.cal
@@ -1,4 +1,4 @@
-SofSci_TOF2INV_V -8.04267 0.00737389
+SofSci_TOF2INV_V -8.0425 0.00737389
SofSci_LENGTH_S2 135.614
-SofSci_DET1_POSPAR 90 54.6
+SofSci_DET1_POSPAR 45.0 80.0
SofSci_DET2_POSPAR 943.287 86.652
diff --git a/Projects/s455/calibration/SofTrim/SofTrim_SectionAlign.cal b/Projects/s455/calibration/SofTrim/SofTrim_SectionAlign.cal
index 292e1b28b..75ecd8429 100644
--- a/Projects/s455/calibration/SofTrim/SofTrim_SectionAlign.cal
+++ b/Projects/s455/calibration/SofTrim/SofTrim_SectionAlign.cal
@@ -1,3 +1,3 @@
-SofTrim_SEC1_ALIGN 755.97 0.987851
+SofTrim_SEC1_ALIGN -348.015 1.03199
SofTrim_SEC2_ALIGN 0 1
-SofTrim_SEC3_ALIGN -1674.22 1.15275
+SofTrim_SEC3_ALIGN 2496.29 0.984162
diff --git a/Projects/s455/configs/ConfigSofSci.dat b/Projects/s455/configs/ConfigSofSci.dat
index 0ac36bed1..fb8753a90 100644
--- a/Projects/s455/configs/ConfigSofSci.dat
+++ b/Projects/s455/configs/ConfigSofSci.dat
@@ -3,6 +3,6 @@ ConfigSofSci
DET1_POSNS_MAX 20
DET2_POSNS_MIN -20
DET2_POSNS_MAX 20
- RAWTOF_MIN 1600
- RAWTOF_MAX 1700
+ RAWTOF_MIN 1628
+ RAWTOF_MAX 1636
diff --git a/Projects/s455/macro/GetChargeDist.C b/Projects/s455/macro/GetChargeDist.C
new file mode 100644
index 000000000..b81f469e0
--- /dev/null
+++ b/Projects/s455/macro/GetChargeDist.C
@@ -0,0 +1,352 @@
+TChain* chain=NULL;
+
+TCutG* cut_Pb1[13];
+TCutG* cut_Pb2[13];
+TCutG* cut_C[13];
+
+int A;
+int Z;
+NPL::Particle* iso;
+
+TH2F* hbid;
+TH1F* hsumpb;
+TH1F* hsumc;
+TH1F* hzpb;
+TH1F* hzc;
+TH1F* hzcscaled;
+TH1F* hsumcscaled;
+TH1F* hzsub;
+
+void LoadRootFiles(string iso);
+
+
+//////////////////////////////////////
+void LoadCuts(){
+
+ TString element[13] = {"180Hg_1", "180Hg_2", "182Hg", "184Hg", "187Pb", "189Pb", "175Pt", "204Fr", "207Fr", "199At_1", "199At_2", "197At", "216Th"};
+
+
+ TFile* file;
+ TString filename1[13];
+ TString filename2[13];
+ TString filename3[13];
+ for(int i=0; i<13; i++){
+ filename1[i]= "cuts/"+element[i]+"/cut_Pb1.root";
+ filename2[i]= "cuts/"+element[i]+"/cut_Pb2.root";
+ filename3[i]= "cuts/"+element[i]+"/cut_C.root";
+
+ file = new TFile(filename1[i],"read");
+ cut_Pb1[i] = (TCutG*) file->FindObjectAny("cut_Pb1");
+ cut_Pb1[i]->SetName(element[i]+"_Pb1");
+
+ file = new TFile(filename2[i],"read");
+ cut_Pb2[i] = (TCutG*) file->FindObjectAny("cut_Pb2");
+ cut_Pb2[i]->SetName(element[i]+"_Pb2");
+
+ file = new TFile(filename3[i],"read");
+ cut_C[i] = (TCutG*) file->FindObjectAny("cut_C");
+ cut_C[i]->SetName(element[i]+"_C");
+ }
+}
+
+//////////////////////////////////////
+void GetChargeDist(string nucleus_name){
+
+ iso = new NPL::Particle(nucleus_name);
+ A = iso->GetA();
+ Z = iso->GetZ();
+ double AoQ=(double)A/(double)Z;
+
+ cout << "*******************************" << endl;
+ cout << "Nucleus name = " << nucleus_name << endl;
+ cout << "A= " << A << endl;
+ cout << "Z= " << Z << endl;
+ cout << "AoQ= " << AoQ << endl;
+
+ LoadRootFiles(nucleus_name);
+ LoadCuts();
+
+
+ //*** opening output root file ***//
+ TString output_filename = "root/charge_dist_" + nucleus_name + ".root";
+ TFile* ofile = new TFile(output_filename, "recreate");
+ cout << "Creating file: " << output_filename << endl;
+ ofile->cd();
+
+ hbid = new TH2F("hbid","hbid",1000,2.2,2.5,1000,70,95);
+ hsumpb = new TH1F("hsumpb","hsumpb",1000,50,95);
+ hsumc = new TH1F("hsumc","hsumc",1000,50,95);
+ hzpb = new TH1F("hzpb","hzpb",1000,20,60);
+ hzc = new TH1F("hzc","hzc",1000,20,60);
+ hzcscaled = new TH1F("hzcscaled","hzcscaled",1000,20,60);
+ hsumcscaled = new TH1F("hsumcscaled","hsumcscaled",1000,50,95);
+ hzsub = new TH1F("hzsub","hzsub",1000,20,60);
+
+
+ //*** Defining beam cut ***//
+ TCutG* cut_beam = new TCutG("cut_beam");
+ cut_beam->SetName("cut_beam");
+ cut_beam->SetVarX("fBeam_AoQ");
+ cut_beam->SetVarY("fBeam_Z");
+ double Rx=0.005;
+ double Ry=0.5;
+ AoQ = AoQ;// + 0.001;
+ double Zcenter = Z;// - 0.5;
+ for(int i=0; i<360; i++){
+ double x = AoQ+Rx*cos(i*TMath::Pi()/180);
+ double y = Zcenter+Ry*sin(i*TMath::Pi()/180);
+ cut_beam->SetPoint(i,x,y);
+ }
+ cut_beam->SetLineWidth(2);
+ cut_beam->SetLineColor(2);
+
+
+
+ int nentries = chain->GetEntries();
+
+ //*** Init input branch ***//
+ int RunID;
+ TSofBeamID* SofBeam = new TSofBeamID();
+ TSofAtPhysics* SofAt = new TSofAtPhysics();
+ TSofFissionFragment* SofFF = new TSofFissionFragment();
+
+ chain->SetBranchAddress("RunID",&RunID);
+ chain->SetBranchAddress("SofBeamID",&SofBeam);
+ chain->SetBranchAddress("SofAt",&SofAt);
+ chain->SetBranchAddress("SofFissionFragment",&SofFF);
+
+ for(int i=0; i<nentries; i++){
+ chain->GetEntry(i);
+
+ if(i%100000==0){
+ cout << "\033[34m\r Processing tree... " << (double)i/nentries*100 << "\% done" << flush;
+ }
+
+ if(SofAt->Energy.size()==3 || SofAt->Energy.size()==4){
+
+ double Zbeam = SofBeam->GetZbeam();
+ double AoQ = SofBeam->GetAoQ();
+
+ double Anode1 = SofAt->Energy[0];
+ double Anode2 = SofAt->Energy[1];
+ double Anode3 = SofAt->Energy[2];
+
+ int k=RunID-1;
+ if(cut_Pb1[k]->IsInside(Anode1,Anode2) || cut_Pb2[k]->IsInside(Anode2,Anode3)){
+ hbid->Fill(AoQ,Zbeam);
+ }
+
+ if(cut_beam->IsInside(AoQ,Zbeam)){
+ //*** Fission Fragment information ***//
+ double Zsum = SofFF->GetZsum();
+ int iZsum = SofFF->GetIntZsum();
+ double Z1 = -1;
+ double Z2 = -1;
+ if(SofFF->GetMult()>0){
+ Z1 = SofFF->GetZ(0);
+ Z2 = SofFF->GetZ(1);
+ }
+
+ //*** Taking event from fission in Lead cathod ***//
+ if(cut_Pb1[k]->IsInside(Anode1,Anode2) || cut_Pb2[k]->IsInside(Anode2,Anode3)){
+ hsumpb->Fill(Zsum);
+ if(iZsum==Z){
+ hzpb->Fill(Z1);
+ hzpb->Fill(Z2);
+
+ hzsub->Fill(Z1);
+ hzsub->Fill(Z2);
+ }
+ }
+
+ //*** Taking event from fission in Carbon cathod ***//
+ if(cut_C[k]->IsInside(Anode2,Anode3)){
+ hsumc->Fill(Zsum);
+ hsumcscaled->Fill(Zsum);
+ if(iZsum==Z){
+ hzc->Fill(Z1);
+ hzc->Fill(Z2);
+
+ hzcscaled->Fill(Z1);
+ hzcscaled->Fill(Z2);
+ }
+ }
+ }
+ }
+ } // end of loop ovent nentries
+
+ double integral_pb = hsumpb->Integral(hsumpb->GetXaxis()->FindBin(65.5), hsumpb->GetXaxis()->FindBin(Z-1.5));
+ double integral_c = hsumc->Integral(hsumc->GetXaxis()->FindBin(65.5), hsumc->GetXaxis()->FindBin(Z-1.5));
+ double ratio = integral_pb/integral_c;
+ if(ratio>0){
+ cout << "*** ratio= " << ratio << endl;
+ hsumcscaled->Scale(ratio);
+ hzcscaled->Scale(ratio);
+
+ hzsub->Add(hzcscaled,-1);
+ }
+ cout << endl;
+ hbid->Write();
+ hsumpb->Write();
+ hsumc->Write();
+ hzpb->Write();
+ hzc->Write();
+ hzcscaled->Write();
+ hsumcscaled->Write();
+ hzsub->Write();
+ cut_beam->Write();
+}
+
+
+//////////////////////////////////////
+void LoadRootFiles(string iso){
+ chain = new TChain("PhysicsTree");
+
+ //*** Th isotopes ***//
+ if(iso=="217Th" || iso=="216Th" || iso=="209Ra" || iso=="210Ra" || iso=="211Ra" || iso=="212Ra" || iso=="213Ra" || iso=="214Ra"){
+ chain->Add("../../../Outputs/Analysis/216Th_1.root");
+ }
+
+ //*** Fr isotopes ***//
+ if(iso=="203Fr" || iso=="204Fr" || iso=="205Fr" || iso=="206Fr"){
+ chain->Add("../../../Outputs/Analysis/204Fr.root");
+ }
+ if(iso=="205Fr" || iso=="206Fr" || iso=="207Fr" || iso=="208Fr"){
+ chain->Add("../../../Outputs/Analysis/207Fr_2.root");
+ }
+ if(iso=="208Fr" || iso=="209Fr" || iso=="210Fr" || iso=="211Fr"){
+ chain->Add("../../../Outputs/Analysis/216Th_1.root");
+ }
+
+ //*** Rn isotopes ***//
+ if(iso=="200Rn" || iso=="201Rn" || iso=="202Rn" || iso=="203Rn" || iso=="204Rn" || iso=="205Rn"){
+ chain->Add("../../../Outputs/Analysis/204Fr.root");
+ }
+ if(iso=="203Rn" || iso=="204Rn" || iso=="205Rn" || iso=="206Rn"){
+ chain->Add("../../../Outputs/Analysis/207Fr.root");
+ }
+
+ //*** At isotopes ***//
+ if(iso=="197At" || iso=="198At" || iso=="199At" || iso=="200At" || iso=="201At" || iso=="202At"){
+ chain->Add("../../../Outputs/Analysis/204Fr.root");
+ }
+ if(iso=="200At" || iso=="201At" || iso=="202At" || iso=="203At"){
+ chain->Add("../../../Outputs/Analysis/207Fr.root");
+ }
+ if(iso=="198At" || iso=="199At" || iso=="200At" || iso=="201At"){
+ chain->Add("../../../Outputs/Analysis/199At_1.root");
+ chain->Add("../../../Outputs/Analysis/199At_2.root");
+ }
+ if(iso=="196At" || iso=="197At" || iso=="198At" || iso=="199At"){
+ chain->Add("../../../Outputs/Analysis/197At.root");
+ }
+
+ //*** Po isotopes ***//
+ if(iso=="192Po" || iso=="193Po" || iso=="194Po" || iso=="195Po" || iso=="196Po" || iso=="197Po"){
+ chain->Add("../../../Outputs/Analysis/197At.root");
+ }
+ if(iso=="195Po" || iso=="196Po" || iso=="197Po" || iso=="198Po" || iso=="199Po" || iso=="200Po"){
+ chain->Add("../../../Outputs/Analysis/199At_1.root");
+ chain->Add("../../../Outputs/Analysis/199At_2.root");
+ }
+
+ //*** Bi isotopes ***//
+ if(iso=="189Bi" || iso=="190Bi" || iso=="191Bi" || iso=="192Bi" || iso=="193Bi" || iso=="194Bi" || iso=="195Bi"){
+ chain->Add("../../../Outputs/Analysis/197At.root");
+ }
+ if(iso=="192Bi" || iso=="193Bi" || iso=="194Bi" || iso=="195Bi" || iso=="196Bi" || iso=="197Bi"){
+ chain->Add("../../../Outputs/Analysis/199At_1.root");
+ chain->Add("../../../Outputs/Analysis/199At_2.root");
+ }
+
+ //*** Pb isotopes ***//
+ if(iso=="186Pb" || iso=="187Pb" || iso=="188Pb" || iso=="189Pb" || iso=="190Pb" || iso=="191Pb" || iso=="192Pb"){
+ chain->Add("../../../Outputs/Analysis/197At.root");
+ }
+ if(iso=="191Pb" || iso=="192Pb" || iso=="193Pb" || iso=="194Pb"){
+ chain->Add("../../../Outputs/Analysis/199At_1.root");
+ chain->Add("../../../Outputs/Analysis/199At_2.root");
+ }
+ if(iso=="188Pb" || iso=="189Pb" || iso=="190Pb" || iso=="191Pb"){
+ chain->Add("../../../Outputs/Analysis/189Pb.root");
+ }
+ if(iso=="186Pb" || iso=="187Pb" || iso=="188Pb" || iso=="189Pb"){
+ chain->Add("../../../Outputs/Analysis/187Pb.root");
+ }
+
+ //*** Tl isotopes ***//
+ if(iso=="186Tl" || iso=="187Tl" || iso=="188Tl" || iso=="189Tl" || iso=="190Tl"){
+ chain->Add("../../../Outputs/Analysis/197At.root");
+ }
+ if(iso=="185Tl" || iso=="186Tl" || iso=="187Tl" || iso=="188Tl" || iso=="189Tl"){
+ chain->Add("../../../Outputs/Analysis/189Pb.root");
+ }
+ if(iso=="183Tl" || iso=="184Tl" || iso=="185Tl" || iso=="186Tl" || iso=="187Tl"){
+ chain->Add("../../../Outputs/Analysis/187Pb.root");
+ }
+ if(iso=="185Tl" || iso=="186Tl" || iso=="187Tl"){
+ chain->Add("../../../Outputs/Analysis/184Hg.root");
+ }
+ if(iso=="183Tl" || iso=="184Tl" || iso=="185Tl"){
+ chain->Add("../../../Outputs/Analysis/182Hg.root");
+ }
+
+ //*** Hg isotopes ***//
+ if(iso=="182Hg" || iso=="183Hg" || iso=="184Hg" || iso=="185Hg" || iso=="186Hg" || iso=="187Hg"){
+ chain->Add("../../../Outputs/Analysis/189Pb.root");
+ }
+ if(iso=="180Hg" || iso=="181Hg" || iso=="182Hg" || iso=="183Hg" || iso=="184Hg" || iso=="185Hg"){
+ chain->Add("../../../Outputs/Analysis/187Pb.root");
+ }
+ if(iso=="183Hg" || iso=="184Hg" || iso=="185Hg" || iso=="186Hg"){
+ chain->Add("../../../Outputs/Analysis/184Hg.root");
+ }
+ if(iso=="181Hg" || iso=="182Hg" || iso=="183Hg" || iso=="184Hg"){
+ chain->Add("../../../Outputs/Analysis/182Hg.root");
+ }
+ if(iso=="178Hg" || iso=="179Hg" || iso=="180Hg" || iso=="181Hg"){
+ chain->Add("../../../Outputs/Analysis/180Hg_1.root");
+ chain->Add("../../../Outputs/Analysis/180Hg_2.root");
+ }
+
+ //*** Au isotopes ***//
+ if(iso=="177Au" || iso=="178Au" || iso=="179Au"){
+ chain->Add("../../../Outputs/Analysis/175Pt.root");
+ }
+ if(iso=="180Au" || iso=="181Au" || iso=="182Au"){
+ chain->Add("../../../Outputs/Analysis/189Pb.root");
+ }
+ if(iso=="179Au" || iso=="180Au" || iso=="181Au" || iso=="182Au"){
+ chain->Add("../../../Outputs/Analysis/187Pb.root");
+ }
+ if(iso=="180Au" || iso=="181Au" || iso=="182Au"){
+ chain->Add("../../../Outputs/Analysis/184Hg.root");
+ }
+ if(iso=="178Au" || iso=="179Au" || iso=="180Au" || iso=="181Au" || iso=="182Au"){
+ chain->Add("../../../Outputs/Analysis/182Hg.root");
+ }
+ if(iso=="177Au" || iso=="178Au" || iso=="179Au" || iso=="180Au"){
+ chain->Add("../../../Outputs/Analysis/180Hg_1.root");
+ chain->Add("../../../Outputs/Analysis/180Hg_2.root");
+ }
+
+ //*** Pt isotopes ***//
+ if(iso=="173Pt" || iso=="174Pt" || iso=="175Pt" || iso=="176Pt" || iso=="177Pt" || iso=="178Pt"){
+ chain->Add("../../../Outputs/Analysis/175Pt.root");
+ }
+ if(iso=="174Pt" || iso=="175Pt" || iso=="176Pt" || iso=="177Pt" || iso=="178Pt"){
+ chain->Add("../../../Outputs/Analysis/180Hg_1.root");
+ chain->Add("../../../Outputs/Analysis/180Hg_2.root");
+ }
+
+ //*** Ir isotopes ***//
+ if(iso=="172Ir" || iso=="173Ir" || iso=="174Ir" || iso=="175Ir" || iso=="176It"){
+ chain->Add("../../../Outputs/Analysis/175Pt.root");
+ chain->Add("../../../Outputs/Analysis/180Hg_1.root");
+ chain->Add("../../../Outputs/Analysis/180Hg_2.root");
+ }
+
+ cout << "Number of entries: " << chain->GetEntries() << endl;
+}
+
+
--
GitLab