rtang/ANASEN_analysis
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

increase smudge

Browse Source
This commit is contained in:
James Szalkie 2026-09-14 13:21:55 -04:00
parent c6c8be8cde
commit e8a43509b9
3 changed files with 76 additions and 23 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

79
Armory/anasenMS.cpp
View File

@ -70,7 +70,7 @@ bool IsDeadCathode(int id){
bool IsDeadSX3(int id){
static std::set<int> dead = {
0, 2, 4, 5, 6, 8, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23
//0, 2, 4, 5, 6, 8, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23
}; // add dead SX3 IDs here, 0-23 1,7,9,3
return dead.count(id);
}
@ -91,7 +91,7 @@ bool IsDeadSX3FrontDnChannel(int sx3ID, int chDn){
bool IsDeadSX3BackChannel(int sx3ID, int chBk){
static std::set<std::pair<int, int>> dead = {
{1, 10}
//{1, 10}
// {sx3ID, back-channel}
};
return dead.count({sx3ID, chBk});
@ -436,6 +436,11 @@ int main(int argc, char **argv){
tree1->Branch("sigma_a", &sigma_a, "sigma_a/D");
tree1->Branch("theta_recon", &theta_recon, "theta_recon/D");
double vX_recon, vY_recon, vZ_recon;
tree1->Branch("vX_recon", &vX_recon, "vX_recon/D");
tree1->Branch("vY_recon", &vY_recon, "vY_recon/D");
tree1->Branch("vZ_recon", &vZ_recon, "vZ_recon/D");
// reconstructed vertex Z from PW fit
double z0;
tree1->Branch("z0", &z0, "reconstucted_Z/D");
@ -483,6 +488,17 @@ int main(int argc, char **argv){
detZ = TMath::QuietNaN();
hitTheta = TMath::QuietNaN();
originalEnergy = TMath::QuietNaN();
aX = TMath::QuietNaN();
aY = TMath::QuietNaN();
aZ = TMath::QuietNaN();
cX = TMath::QuietNaN();
cY = TMath::QuietNaN();
cZ = TMath::QuietNaN();
theta_recon = TMath::QuietNaN();
vX_recon = TMath::QuietNaN();
vY_recon = TMath::QuietNaN();
vZ_recon = TMath::QuietNaN();
Ex_recon = TMath::QuietNaN();
sx3ID = -1;
sx3Up = -1;
sx3Dn = -1;
@ -511,8 +527,8 @@ int main(int argc, char **argv){
double vertexRangeX = std::sqrt((vertexXRange[1] * vertexXRange[1]) + (sigmaX_mm * sigmaX_mm));
double vertexRangeY = std::sqrt((vertexYRange[1] * vertexYRange[1]) + (sigmaX_mm * sigmaX_mm));
//std::cout << "vertexRangeX: " << vertexRangeX << ", vertexRangeY: " << vertexRangeY << std::endl;
vertexX = gRandom->Gaus(0, vertexRangeX); // mean and standard deviation
vertexY = gRandom->Gaus(0, vertexRangeY); // mean and standard deviation
vertexX = gRandom->Gaus(-5, vertexRangeX); // mean and standard deviation
vertexY = gRandom->Gaus(7, vertexRangeY); // mean and standard deviation
//vertexZ = (vertexZRange[1]- vertexZRange[0])*gRandom->Rndm() + vertexZRange[0];
@ -529,7 +545,8 @@ int main(int argc, char **argv){
//KEA = beamEnergy / beamA;
//KEA = gRandom->Uniform(0, beamE);
//sigma_a = std::max(0.0, sigmaABeam->Eval(beamEnergy));
transfer.SetIncidentEnergyAngle(KEA, 0, 0); //arguments are (kinetic energy, polar angle, azimuthal angle) of the incident particle in the lab frame
double beamTheta = gRandom->Gaus(0, 10 * TMath::DegToRad());
transfer.SetIncidentEnergyAngle(KEA, beamTheta, 0); //arguments are (kinetic energy, polar angle, azimuthal angle) of the incident particle in the lab frame
transfer.CalReactionConstant();
// isotropic CM direction
@ -662,14 +679,9 @@ int main(int argc, char **argv){
//TVector3 hitSigma = sx3->GetHitPosWithSigma(sigmaSX3_W, sigmaSX3_L);
TVector3 hitSigma = sx3->GetHitPosWithSigma(0, 0);
double hitX = std::lround(hitSigma.X());
double hitY = std::lround(hitSigma.Y());
double hitZ = std::lround(hitSigma.Z());
originalEnergy = CalculateOriginalEnergy(hitX, hitY, hitZ,
0, 0, std::lround(vertexZ / 10) * 10,
b, "He", Esx3,
distance_sx3);
double hitX = hitSigma.X();
double hitY = hitSigma.Y();
double hitZ = hitSigma.Z();
double distance_A;
Eanode = CalculateEnergyLoss(vertexX, vertexY, vertexZ,
@ -698,12 +710,41 @@ int main(int argc, char **argv){
//tree1->Fill();
Kinematics apkin_27Al(26.981538408,4.00260325413,1.00782503224,29.973770136,beamEnergy/26.981538408); //m3 is proton
//reconstruct theta using anode and sx3 positions, apply a sigma to the anode locations
double aX_sigma = gRandom->Gaus(aX, 3);
double aY_sigma = gRandom->Gaus(aY, 3);
double aZ_sigma = gRandom->Gaus(aZ, 3);
theta_recon = std::asin(TVector3(sx3X - aX_sigma, sx3Y - aY_sigma, 0).Mag() / TVector3(sx3X - aX_sigma, sx3Y - aY_sigma, sx3Z - aZ_sigma).Mag()) * 180.0 / TMath::Pi();
Ex_recon = apkin_27Al.getExc(originalEnergy, theta_recon);
if (!std::isnan(aX) && !std::isnan(aY) && !std::isnan(aZ) && aZ != 54321 && (aX != 0 || aY != 0 || aZ != 0)) {
//reconstruct theta using anode and sx3 positions, apply a sigma to the anode locations
double aX_sigma = gRandom->Gaus(aX, 0);
double aY_sigma = gRandom->Gaus(aY, 0);
double aZ_sigma = gRandom->Gaus(aZ, 0);
theta_recon = TVector3(sx3X - aX_sigma, sx3Y - aY_sigma, sx3Z - aZ_sigma).Theta() * 180.0 / TMath::Pi();
//calulate vertex like in experiment
TVector3 x2f(aX_sigma, aY_sigma, aZ_sigma);
TVector3 x1(sx3X, sx3Y, sx3Z);
TVector3 v = x2f - x1;
double t_minimum = -1.0*(x1.X()*v.X()+x1.Y()*v.Y())/(v.X()*v.X()+v.Y()*v.Y());
TVector3 r_rhoMin_fix = x1 + t_minimum*v;
originalEnergy = CalculateOriginalEnergy(hitX, hitY, hitZ,
r_rhoMin_fix.X(), r_rhoMin_fix.Y(), r_rhoMin_fix.Z(),
//0, 0, r_rhoMin_fix.Z(),
//0, 0, vertexZ,
b, "He", Esx3,
distance_sx3);
vX_recon = r_rhoMin_fix.X();
vY_recon = r_rhoMin_fix.Y();
vZ_recon = r_rhoMin_fix.Z();
Ex_recon = apkin_27Al.getExc(originalEnergy, theta_recon);
/*Checklist: anode smudge, sx3 smudge, beam position off axis, beam angle*/
} else {
theta_recon = TMath::QuietNaN();
vX_recon = TMath::QuietNaN();
vY_recon = TMath::QuietNaN();
vZ_recon = TMath::QuietNaN();
Ex_recon = TMath::QuietNaN();
originalEnergy = TMath::QuietNaN();
}
//EBeam_Kin_gs = apkin_27Al.getEbeam_givenQ(Esx3, 0.0, thetab);
//EBeam_Kin_2_2 = apkin_27Al.getEbeam_givenQ(Esx3, 2.2, thetab);
//EBeam_Kin_3_4 = apkin_27Al.getEbeam_givenQ(Esx3, 3.4, thetab);

2
Armory/anasen_anode_cathode_hyperboloids.h
View File

@ -31,7 +31,7 @@ inline std::tuple<TVector3,TVector3,double> find_PC_PathLength(const TVector3& x
TVector3 anode_intersect = onesheet_hyperboloid_intersect(32.0429,301.895);
TVector3 cathode_intersect = onesheet_hyperboloid_intersect(37.239045,301.895);
TVector3 gw_intersect = onesheet_hyperboloid_intersect(27.712,301.895);
if(anode_intersect.Z()!=54321 && cathode_intersect.Z()!=54321)
if(anode_intersect.Z()!=54321 && cathode_intersect.Z()!=54321 && gw_intersect.Z()!=54321)
return std::tuple(cathode_intersect,gw_intersect,(cathode_intersect-gw_intersect).Mag()*0.1);
else
return std::tuple(TVector3(0,0,0), TVector3(0,0,0), 54321);

18
Armory/testMacro.C
View File

@ -8,6 +8,7 @@
TH2F *h2 = new TH2F("EPCvEsx3", "EPC x sin(thetab) vs Esx3;Esx3;EPC * sin(thetab * TMath::DegToRad())", 200, 0, 0, 200, 0, 0); //arguments are (name, title, nbinsX, xlow, xup, nbinsY, ylow, yup)
tree1->Draw("EPC*sin(thetab * TMath::DegToRad()):Esx3 >> EPCvEsx3", "Esx3 > 0", "colz");
/*
new TCanvas();
TH2F *h3 = new TH2F(
"hExReconVsVZ",
@ -23,7 +24,7 @@
100, 0, 0,
100, -1, 0
);
tree1->Draw("Ex:vZ>>hExVsVZ", "sx3ID >= 0", "box same");
tree1->Draw("Ex:vZ>>hExVsVZ", "sx3ID >= 0", "box same");*/
/*
new TCanvas();
tree1->Draw("originalEnergy:Tb", "sx3ID >=0", "");
@ -31,8 +32,19 @@
f_diagonal->SetLineColor(kRed);
f_diagonal->SetLineStyle(2);
f_diagonal->Draw("same");*/
/*
TH2F *h5 = new TH2F("Excitation Energy vs vZ", "Excitation Energy vs vZ;vZ (mm);Excitation energy (MeV)", 300, 0, 0, 300, 0, 0); //arguments are (name, title, nbinsX, xlow, xup, nbinsY, ylow, yup)
tree1->Draw("Ex_recon:vZ >> hExReconVsVZ", "sx3ID >= 0", "colz");
tree1->Draw("Ex:vZ >> h5", "sx3ID >= 0 && EPC >= 0.03", "same");
*/
new TCanvas("c1", "Canvas", 900, 600);
tree1->Draw("Ex_recon:vZ_recon>>h2(800,-500,300,110,-1,10)", "sx3ID >=0", "colz");
tree1->Draw("Ex:vZ", "sx3ID >=0", "same");
c1->SaveAs("Ex_recon.png");
new TCanvas("c2", "Canvas", 900, 600);
tree1->Draw("theta_recon:thetab", "sx3ID >=0", "");
c2->SaveAs("theta_recon_vs_thetab.png");
TH2F *h5 = new TH2F("Excitation Energy vs vZ", "Excitation Energy vs vZ;vZ (mm);Excitation energy (MeV)", 200, 0, 0, 200, 0, 0); //arguments are (name, title, nbinsX, xlow, xup, nbinsY, ylow, yup)
tree1->Draw("Ex:vZ >> Excitation Energy vs vZ", "sx3ID >= 0", "colz");
}