Skip to content
Snippets Groups Projects
Commit 8abed143 authored by Charlie Paxman's avatar Charlie Paxman
Browse files

* e793s - Autoselect solid angle sim

parent 7b620b9d
No related branches found
No related tags found
No related merge requests found
Pipeline #203868 passed
Stage: build-NPLib
Stage: build-NPSimulation
Stage: test
Hide whitespace changes
Inline Side-by-side
Projects/e793s/macro/CS2.h
+ 39
− 21
View file @ 8abed143
......@@ -21,10 +21,14 @@ bool loud = 1;
/* Scale method toggle */
bool scaleTogether = 1;
/* String for image */
/* Strings for image */
string orbitalname;
string orbital;
/* Strings for SolidAngle input file */
string statename;
string inputdate;
////////////////////////////////////////////////////////////////////////////////
void canclone(TCanvas* major, int padNum, string name){
......@@ -84,22 +88,6 @@ void CS_Diagnosis(){
////////////////////////////////////////////////////////////////////////////////
void CS(){
/* Overload function */
/*
cout << " Inputs:\n Experimental...\n\t - Energy of state\n\t - Spin of state" << endl;
cout << " Theory...\n\t - Orbital l\n\t - Orbital j\n\t - Orbital n\n\n" << endl;
cout << " 0f5/2 | ----- | --?-- | l=3, j=2.5, n=0" << endl;
cout << " | | |" << endl;
cout << " 1p1/2 | ----- | --?-- | l=1, j=0.5, n=1" << endl;
cout << " 1p3/2 | ----- | --?-- | l=1, j=1.5, n=1" << endl;
cout << " 0f7/2 | ----- | ===== |" << endl;
cout << " | | |" << endl;
cout << " 0d3/2 | ----- | ===== |" << endl;
cout << " 1s1/2 | --x-- | ===== |" << endl;
cout << " | | |" << endl;
cout << " | p | n |" << endl;
*/
cout << "- CS(stateE, stateSp, orb_l, orb_j, nodes) "<< endl;
cout << "---- 0.143, p3/2 = CS(0.143, 2, 1, 1.5) "<< endl;
cout << "---- 0.279, p3/2 = CS(0.279, 2, 1, 1.5) "<< endl;
......@@ -113,8 +101,6 @@ void CS(){
cout << "---- 3.8 , f5/2 = CS(3.8 , 3, 3, 2.5) "<< endl;
cout << "---- 4.1 , f5/2 = CS(4.1 , 3, 3, 2.5) "<< endl;
cout << "---- 4.4 , f5/2 = CS(4.4 , 3, 3, 2.5) "<< endl;
}
////////////////////////////////////////////////////////////////////////////////
......@@ -129,7 +115,7 @@ void CS(double Energy, double Spin, double spdf, double angmom){
//double ElasticNorm = 3.7, ElasticNormErr = 0.000; // Estimated 'goal' normalization
double ElasticNorm = 0.000220, ElasticNormErr = 0.000; //14Oct22
inputdate = "18Oct22";
orbitalname.clear();
orbital.clear();
......@@ -195,6 +181,9 @@ void CS(double Energy, double Spin, double spdf, double angmom){
cout << "Identified as state #" << i << ", E = " << means[i] << endl;
indexE = i;
found = 1;
stringstream ss;
ss << setfill('0') << setw(4) << (int)(means[i]*1000.);
statename = ss.str();
}
}
if(!found){
......@@ -210,10 +199,39 @@ void CS(double Energy, double Spin, double spdf, double angmom){
//auto file = new TFile("../SolidAngle_HistFile_19Jul22_47Kdp_1p981.root");
//auto file = new TFile("../SolidAngle_HistFile_19Jul22_47Kdp_4p393.root");
//auto file = new TFile("../SolidAngle_HistFile_30Jul22_47Kdp_0p000_ThetaBin0p5.root");
auto file = new TFile("../SolidAngle_HistFile_10Aug22_TrueStripRemoval.root");
//auto file = new TFile("SolidAngle_HistFiles/SolidAngle_HistFile_10Aug22_TrueStripRemoval.root");
string backupFileName = "SolidAngle_HistFiles/SolidAngle_HistFile_10Aug22_TrueStripRemoval.root";
string saFileName = "SolidAngle_HistFiles/SAHF_";
saFileName.append(inputdate);
saFileName.append("_47Kdp_");
saFileName.append(statename);
saFileName.append(".root");
TFile* file;
ifstream infile(saFileName);
ifstream backup(backupFileName);
if(infile.good()){
cout << "Opening file " << saFileName << endl;
file = TFile::Open(saFileName.c_str());
} else if (backup.good()){
cout << BOLDRED << "FAILED TO OPEN " << saFileName << endl;
cout << "Open standard backup file " << backupFileName << endl;
cout << RESET;
file = TFile::Open(backupFileName.c_str());
} else {
cout << BOLDRED << "FAILED TO OPEN MAIN OR BACKUP SOLID ANGLE FILE" << endl;
cout << RED << "Check SolidAngle file exists..." << RESET << endl;
}
//cout << "MADE IT OUT" << endl;
//auto file = new TFile("../SolidAngle_HistFile_New.root");
/* ADD OPTION TO CHANGE SOLID ANGLE FILE DEPENDING ON PEAK!!!!*/
TH1F* SolidAngle = (TH1F*) file->FindObjectAny("SolidAngle_Lab_MG");
//cout << BLUE << "MADE IT HERE" << endl;
TCanvas* c_SolidAngle = new TCanvas("c_SolidAngle","c_SolidAngle",1000,1000);
SolidAngle->Draw("HIST");
SolidAngle->GetXaxis()->SetRangeUser(100.,160.);
......
Projects/e793s/macro/CS2_dt.h
+ 34
− 48
View file @ 8abed143
......@@ -25,6 +25,9 @@ bool scaleTogether = 1;
string orbitalname;
string orbital;
/* Strings for SolidAngle input file */
string statename;
string inputdate;
////////////////////////////////////////////////////////////////////////////////
void canclone(TCanvas* major, int padNum, string name){
......@@ -84,37 +87,7 @@ void CS_Diagnosis(){
////////////////////////////////////////////////////////////////////////////////
void CS(){
/* Overload function */
/*
cout << " Inputs:\n Experimental...\n\t - Energy of state\n\t - Spin of state" << endl;
cout << " Theory...\n\t - Orbital l\n\t - Orbital j\n\t - Orbital n\n\n" << endl;
cout << " 0f5/2 | ----- | --?-- | l=3, j=2.5, n=0" << endl;
cout << " | | |" << endl;
cout << " 1p1/2 | ----- | --?-- | l=1, j=0.5, n=1" << endl;
cout << " 1p3/2 | ----- | --?-- | l=1, j=1.5, n=1" << endl;
cout << " 0f7/2 | ----- | ===== |" << endl;
cout << " | | |" << endl;
cout << " 0d3/2 | ----- | ===== |" << endl;
cout << " 1s1/2 | --x-- | ===== |" << endl;
cout << " | | |" << endl;
cout << " | p | n |" << endl;
*/
cout << "- CS(stateE, stateSp, orb_l, orb_j, nodes) "<< endl;
cout << "---- 0.143, p3/2 = CS(0.143, 2, 1, 1.5) "<< endl;
cout << "---- 0.279, p3/2 = CS(0.279, 2, 1, 1.5) "<< endl;
cout << "---- 0.728, f7/2 = CS(0.728, 3, 3, 3.5) "<< endl;
cout << "---- 0.968, p1/2 = CS(0.968, 0, 1, 0.5) "<< endl;
cout << "---- 1.410, p3/2 = CS(1.410, 1, 1, 1.5) "<< endl;
cout << "---- 1.981, p3/2 = CS(1.981, 1, 1, 0.5) "<< endl;
cout << "---- 2.410, p3/2 = CS(2.410, 0, 1, 0.5) "<< endl;
cout << "---- 3.2 , f7/2 = CS(3.2 , 3, 3, 3.5) "<< endl;
cout << "---- 3.6 , f5/2 = CS(3.6 , 3, 3, 2.5) "<< endl;
cout << "---- 3.8 , f5/2 = CS(3.8 , 3, 3, 2.5) "<< endl;
cout << "---- 4.1 , f5/2 = CS(4.1 , 3, 3, 2.5) "<< endl;
cout << "---- 4.4 , f5/2 = CS(4.4 , 3, 3, 2.5) "<< endl;
}
////////////////////////////////////////////////////////////////////////////////
......@@ -122,15 +95,11 @@ void CS(double Energy, double Spin, double spdf, double angmom){
// p3/2 -> spdf = 1, angmom = 1.5
// J0 is incident spin, which is 47K g.s. therefore J0 = 1/2
double J0 = 0.5;
//double ElasticNorm = 5.8, ElasticNormErr = 1.3; // DeuteronNorm in elastics, 5.8 +- 1.3
//double ElasticNorm = 2.70, ElasticNormErr = 0.46; // For thicker
//double ElasticNorm = 2.8, ElasticNormErr = 0.7; // DeuteronNorm in elastics, 2.8 +- 0.7
//double ElasticNorm = 2.363, ElasticNormErr = 0.000; // DeuteronNorm in elastics, reconstructed from determined thickness
//double ElasticNorm = 3.7, ElasticNormErr = 0.000; // Estimated 'goal' normalization
double ElasticNorm = 0.000220, ElasticNormErr = 0.000; //14Oct22
inputdate = "18Oct22";
// Extract spin orbit name
orbitalname.clear();
orbital.clear();
if(spdf==0){
......@@ -165,15 +134,8 @@ void CS(double Energy, double Spin, double spdf, double angmom){
orbital="???";
}
// EDITED UP TO HERE!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
/* Reduce by factor of 10,000 */
//ElasticNorm /= 10000.;
////ElasticNorm /= 1000.;
//ElasticNormErr /= 10000.;
////ElasticNormErr /= 1000.;
// Number of nodes
double nodes;
if(spdf==0){
nodes=1;
}
......@@ -202,6 +164,9 @@ void CS(double Energy, double Spin, double spdf, double angmom){
cout << "Identified as state #" << i << ", E = " << means_dt[i] << endl;
indexE = i;
found = 1;
stringstream ss;
ss << setfill('0') << setw(4) << (int)(means[i]*1000.);
statename = ss.str();
}
}
if(!found){
......@@ -211,9 +176,30 @@ void CS(double Energy, double Spin, double spdf, double angmom){
}
/* Solid Angle (from simulation) */
/* ADD OPTION TO CHANGE SOLID ANGLE FILE DEPENDING ON PEAK!!!!*/
auto file = new TFile("../SolidAngle_HistFile_18Oct22_47Kdt.root");
/* ADD OPTION TO CHANGE SOLID ANGLE FILE DEPENDING ON PEAK!!!!*/
//auto file = new TFile("../SolidAngle_HistFile_18Oct22_47Kdt.root");
string backupFileName = "SolidAngle_HistFiles/SolidAngle_HistFile_18Oct22_47Kdt.root";
string saFileName = "SolidAngle_HistFiles/SAHF_";
saFileName.append(inputdate);
saFileName.append("_47Kdt_");
saFileName.append(statename);
saFileName.append(".root");
TFile* file;
ifstream infile(saFileName);
ifstream backup(backupFileName);
if(infile.good()){
cout << "Opening file " << saFileName << endl;
file = TFile::Open(saFileName.c_str());
} else if (backup.good()){
cout << BOLDRED << "FAILED TO OPEN " << saFileName << endl;
cout << "Open standard backup file " << backupFileName << endl;
cout << RESET;
file = TFile::Open(backupFileName.c_str());
} else {
cout << BOLDRED << "FAILED TO OPEN MAIN OR BACKUP SOLID ANGLE FILE" << endl;
cout << RED << "Check SolidAngle file exists..." << RESET << endl;
}
TH1F* SolidAngle = (TH1F*) file->FindObjectAny("SolidAngle_CM_MM");
TCanvas* c_SolidAngle = new TCanvas("c_SolidAngle","c_SolidAngle",1000,1000);
SolidAngle->Draw("HIST");
......
Projects/e793s/macro/DrawPlots.h
+ 20
− 11
View file @ 8abed143
......@@ -115,8 +115,15 @@ void LoadChain47Kpp(){
//files.push_back("../../../Outputs/Analysis/47Kpp_08Nov_PartI.root");
//files.push_back("../../../Outputs/Analysis/47Kpp_08Nov_PartII.root");
files.push_back("../../../Outputs/Analysis/47Kpp_11Jul22_PartI.root");
files.push_back("../../../Outputs/Analysis/47Kpp_11Jul22_PartII.root");
//files.push_back("../../../Outputs/Analysis/47Kpp_11Jul22_PartI.root");
//files.push_back("../../../Outputs/Analysis/47Kpp_11Jul22_PartII.root");
//files.push_back("../../../Outputs/Analysis/24Oct22_47Kpp_PartI.root");
//files.push_back("../../../Outputs/Analysis/24Oct22_47Kpp_PartII.root");
files.push_back("../../../Outputs/Analysis/25Oct22_47Kpp_PartI.root");
files.push_back("../../../Outputs/Analysis/25Oct22_47Kpp_PartII.root");
chain = Chain("PhysicsTree",files,true);
}
......@@ -384,18 +391,20 @@ void Draw_2DParticleGamma(){
ExEg->GetXaxis()->SetTitle("E_{#gamma} [MeV]");
ExEg->GetYaxis()->SetRangeUser(-1.0,7.0);
ExEg->Draw();
if(reactionName=="47K(d,p)"){
TLine *Sn = new TLine(0,4.644,4.644,4.644);
Sn->SetLineColor(kRed);Sn->SetLineStyle(2);
Sn->Draw("same");
TLatex *TSn = new TLatex(.5,.5,"S_{n}");
TSn->SetTextColor(kRed);
TSn->SetTextSize(0.05);
TSn->SetX(2.50);
TSn->SetY(4.90);
TSn->Draw("same");
}
TLine *XeqY = new TLine(0,0,5,5);
XeqY->SetLineColor(kRed);
XeqY->Draw("same");
TLine *Sn = new TLine(0,4.644,4.644,4.644);
Sn->SetLineColor(kRed);Sn->SetLineStyle(2);
Sn->Draw("same");
TLatex *TSn = new TLatex(.5,.5,"S_{n}");
TSn->SetTextColor(kRed);
TSn->SetTextSize(0.05);
TSn->SetX(2.50);
TSn->SetY(4.90);
TSn->Draw("same");
TLatex *Texeg = new TLatex(.5,.5,"Ex = E_{#gamma}");
Texeg->SetTextColor(kRed);
Texeg->SetTextSize(0.05);
......
Projects/e793s/macro/GausFit.h
+ 2
− 0
View file @ 8abed143
......@@ -94,6 +94,8 @@ vector<double> DoubleGausNumbs(TH1F* hist, double minFit, double maxFit, double
//g2->SetParLimits(2, 0.05, 0.20);
g2->SetParLimits(2, 0.4, 1.0);//FOR ELASTICS
bg->SetParameter(0, 5.0);//FOR ELASTICS
bg->SetParLimits(0, 0.0, 20.0);//FOR ELASTICS
hist->Fit(g1, "WWR", "", minFit, mean1+5);//maxFit);
hist->Fit(g2, "WWR", "", mean2-5, maxFit);//minFit, maxFit);
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment