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

modified: TrackRecon.C the code was NOT dead

Browse Source
This commit is contained in:
Vignesh Sitaraman 2026-09-18 17:48:06 -04:00
parent 397c410320
commit f30e244d27
1 changed files with 45 additions and 17 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

62
TrackRecon.C
View File

@ -63,6 +63,9 @@ double source_vertex = 53.0,
z_entrance = -174.3 - 9.7 - 270.0, z_entrance = -174.3 - 9.7 - 270.0,
dither_sigma = 8.0, dither_sigma = 8.0,
cathode_gain = 1.0, cathode_gain = 1.0,
// Low-band split controls, read by a1c1_solve() in Armory/PCZRecon.h.
// cfrac_split and lowband_rfactor are set per-dataset in Begin();
// missing_fmax is the dead-neighbour cell f-limit in solve_cell().
a1c1_cfrac_split = 0.0, a1c1_cfrac_split = 0.0,
a1c1_missing_fmax = 2.0, a1c1_missing_fmax = 2.0,
a1c1_lowband_rfactor = 0.0, a1c1_lowband_rfactor = 0.0,
@ -79,7 +82,7 @@ double source_vertex = 53.0,
alpha_source_mev = 5.486; alpha_source_mev = 5.486;
// --- Immutable Constants --- // --- Immutable Constants ---
const double qqq_z = 100.0, const double qqq_z = 105.0,
sx3_phi_pitch = 6.5 * (M_PI / 180.0), sx3_phi_pitch = 6.5 * (M_PI / 180.0),
qqq_wedge_pitch = (87.0 / 16.0) * (M_PI / 180.0), qqq_wedge_pitch = (87.0 / 16.0) * (M_PI / 180.0),
qqq_ring_pitch = 48.0 / 16.0; qqq_ring_pitch = 48.0 / 16.0;
@ -271,25 +274,26 @@ struct AAEjectileMasses
const double a1c1_zg[8] = {147.998, 101.946, 59.7634, 19.6965, -19.6965, -59.7634, -101.946, -147.998}; const double a1c1_zg[8] = {147.998, 101.946, 59.7634, 19.6965, -19.6965, -59.7634, -101.946, -147.998};
static const double a1c1_cfmin_17F[7] = {0.20, 0.20, 0.20, 0.20, 0.20, 0.20, 0.20}; static const double a1c1_cfmin_17F[7] = {0.410, 0.420, 0.400, 0.410, 0.410, 0.420, 0.410};
static const double a1c1_k_17F[7] = {0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25}; static const double a1c1_k_17F[7] = {0.060, 0.060, 0.052, 0.059, 0.065, 0.060, 0.060};
static const double a1c1_cfmin_27Al[7] = {0.42, 0.42, 0.42, 0.40, 0.42, 0.43, 0.43}; static const double a1c1_cfmin_27Al[7] = {0.410, 0.420, 0.400, 0.420, 0.400, 0.420, 0.410};
static const double a1c1_k_27Al[7] = {0.06, 0.06, 0.06, 0.06, 0.06, 0.06, 0.06}; static const double a1c1_k_27Al[7] = {0.075, 0.099, 0.066, 0.062, 0.067, 0.099, 0.075};
// low band for 17F data // low band for 17F data
static const double a1c1_cfmin2_17F[7] = {0.10, 0.10, 0.10, 0.10, 0.10, 0.10, 0.10}; static const double a1c1_cfmin2_17F[7] = {0.10, 0.10, 0.10, 0.10, 0.10, 0.10, 0.10};
static const double a1c1_k2_17F[7] = {0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05}; static const double a1c1_k2_17F[7] = {0.05, 0.05, 0.05, 0.05, 0.05, 0.05, 0.05};
static const double a1c1_cfmin2_27Al[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; // no low band static const double a1c1_cfmin2_27Al[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; // no low band
static const double a1c1_k2_27Al[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; static const double a1c1_k2_27Al[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
// active per-cell set, populated by dataset in Begin()
double a1c1_cfmin_cell[7] = {0.410, 0.420, 0.400, 0.410, 0.410, 0.420, 0.410};
double a1c1_k_cell[7] = {0.060, 0.060, 0.052, 0.059, 0.065, 0.060, 0.060};
// Low-band split cells, read by a1c1_solve() in Armory/PCZRecon.h; populated
// per-dataset in Begin() from the a1c1_cfmin2_*/a1c1_k2_* sets above.
double a1c1_cfmin2_cell[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; double a1c1_cfmin2_cell[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
double a1c1_k2_cell[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}; double a1c1_k2_cell[7] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
// active per-cell set, populated by dataset in Begin()
double a1c1_cfmin_cell[7] = {0.20, 0.20, 0.20, 0.20, 0.20, 0.20, 0.20};
double a1c1_k_cell[7] = {0.25, 0.25, 0.25, 0.25, 0.25, 0.25, 0.25};
static std::vector<int> a1c1_dead_anode_17F = {9, 12}; // 1 can be recovered static std::vector<int> a1c1_dead_anode_17F = {9, 12}; // 1 can be recovered
static std::vector<int> a1c1_dead_cathode_17F = {}; // 0,13,15 can be recovered static std::vector<int> a1c1_dead_cathode_17F = {}; // 0,13,15 can be recovered
static std::vector<int> a1c1_dead_anode_27Al = {0, 12, 19}; static std::vector<int> a1c1_dead_anode_27Al = {0, 12, 19};
@ -2005,7 +2009,17 @@ Bool_t TrackRecon::Process(Long64_t entry)
// siPcCoincident(): this used to be the one Si-PC match in the file with a phi // siPcCoincident(): this used to be the one Si-PC match in the file with a phi
// window but no time gate, so "withPC" included phi-aligned but time-accidental // window but no time gate, so "withPC" included phi-aligned but time-accidental
// pairs. // pairs.
auto hasPCCoincidence = [&](const Event &sievent, bool isQQQ) // Per-track PC topology class of the first phi/time-coincident PC hit. The
// anode multiplicity (multi1) is split out because a 2-wire anode cluster
// interpolates between wires (charge-weighted pseudowire) and so resolves z
// better than the wire-pitch-quantized single wire:
// 0 = no PC coincidence
// 1 = A1C0 2 = A2C0 (multi2==0, anode-only)
// 3 = A1C1 4 = A2C1 (multi2==1)
// 5 = A1C2 6 = A2C2 (multi2==2)
// 7 = other (multi2>=3)
// (class>0 reproduces the old hasPCCoincidence bool exactly.)
auto pcTopoClass = [&](const Event &sievent, bool isQQQ) -> int
{ {
double phi_win = isQQQ ? TMath::Pi() / 4.0 : TMath::Pi() / 3.0; double phi_win = isQQQ ? TMath::Pi() / 4.0 : TMath::Pi() / 3.0;
for (const auto &pcevent : PC_Events) for (const auto &pcevent : PC_Events)
@ -2015,21 +2029,35 @@ Bool_t TrackRecon::Process(Long64_t entry)
if (!siPcCoincident(sievent.Time1, pcevent.Time1)) if (!siPcCoincident(sievent.Time1, pcevent.Time1))
continue; continue;
if (TMath::Abs(sievent.pos.DeltaPhi(pcevent.pos)) <= phi_win) if (TMath::Abs(sievent.pos.DeltaPhi(pcevent.pos)) <= phi_win)
return true; {
int c = pcevent.multi2;
if (c > 3)
return 7;
return 1 + 2 * c + (pcevent.multi1 >= 2 ? 1 : 0); // A1C0..A2C2
}
} }
return false; return 0;
}; };
for (const auto &qqqevent : QQQ_Events) for (const auto &qqqevent : QQQ_Events)
{ {
plotter->Fill1D("siE_qqq_calibrated_all", 800, 0, 15, qqqevent.Energy1, "siE"); plotter->Fill1D("siE_qqq_calibrated_all", 800, 0, 15, qqqevent.Energy1, "siE");
bool coinc = hasPCCoincidence(qqqevent, true); int topo = pcTopoClass(qqqevent, true);
plotter->Fill1D(coinc ? "siE_qqq_calibrated_withPC" : "siE_qqq_calibrated_noPC", 800, 0, 15, qqqevent.Energy1, "siE"); plotter->Fill1D(topo > 0 ? "siE_qqq_calibrated_withPC" : "siE_qqq_calibrated_noPC", 800, 0, 15, qqqevent.Energy1, "siE");
// PC acceptance / z-quality mix vs track angle (species compared across files):
// class 0=noPC 1=A1C0 2=A1C1 3=A1C2 4=multi2>=3. theta rel. to source_vertex
// (exact for the point source; a lab-angle proxy for reaction protons).
double th = (qqqevent.pos - beamAxisPoint(source_vertex)).Theta() * 180.0 / M_PI;
plotter->Fill1D("PCtopo_class_qqq", 8, 0, 8, topo + 0.5, "PCtopo");
plotter->Fill2D("PCtopo_class_vs_theta_qqq", 180, 0, 180, 8, 0, 8, th, topo + 0.5, "PCtopo");
} }
for (const auto &sx3event : SX3_Events) for (const auto &sx3event : SX3_Events)
{ {
plotter->Fill1D("siE_sx3_calibrated_all", 800, 0, 15, sx3event.Energy1, "siE"); plotter->Fill1D("siE_sx3_calibrated_all", 800, 0, 15, sx3event.Energy1, "siE");
bool coinc = hasPCCoincidence(sx3event, false); int topo = pcTopoClass(sx3event, false);
plotter->Fill1D(coinc ? "siE_sx3_calibrated_withPC" : "siE_sx3_calibrated_noPC", 800, 0, 15, sx3event.Energy1, "siE"); plotter->Fill1D(topo > 0 ? "siE_sx3_calibrated_withPC" : "siE_sx3_calibrated_noPC", 800, 0, 15, sx3event.Energy1, "siE");
double th = (sx3event.pos - beamAxisPoint(source_vertex)).Theta() * 180.0 / M_PI;
plotter->Fill1D("PCtopo_class_sx3", 8, 0, 8, topo + 0.5, "PCtopo");
plotter->Fill2D("PCtopo_class_vs_theta_sx3", 180, 0, 180, 8, 0, 8, th, topo + 0.5, "PCtopo");
} }
if (doMiscHistograms && ta_foil_run) if (doMiscHistograms && ta_foil_run)