From 96ccdd6e96a64e68f3863010ab378fdf6059b8c4 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Mon, 30 Mar 2026 22:04:51 -0400 Subject: [PATCH 01/10] First attempt at fixing bug --- AtTools/AtELossBetheBloch.cxx | 45 +++++++++++++++---- AtTools/AtELossBetheBloch.h | 4 +- AtTools/AtELossBetheBlochTest.cxx | 74 ++++++++++++++++++++++++++++--- 3 files changed, 108 insertions(+), 15 deletions(-) diff --git a/AtTools/AtELossBetheBloch.cxx b/AtTools/AtELossBetheBloch.cxx index 91c7e283a..90437c094 100644 --- a/AtTools/AtELossBetheBloch.cxx +++ b/AtTools/AtELossBetheBloch.cxx @@ -78,6 +78,24 @@ void AtELossBetheBloch::BuildSpline(double E_min_MeV, double E_max_MeV, int nPoi dXdE[i] = (dedxValues[i] > 0) ? 1.0 / dedxValues[i] : 1e10; fdXdE = tk::spline(energies, dXdE); + + // Build Ω²(E) spline: Bohr range variance accumulated from 0 to E. + // ω²_unit [MeV²/mm] = K · z² · (Z/A) · ρ[g/mm³] · mₑc² + // (kK is in MeV cm²/mol; ×100 → MeV mm²/mol; density g/cm³ → g/mm³ = /1000; combined: /10) + double z = IsElectron() ? 1.0 : fPart_q; + double omega2_unit = kK * z * z * (static_cast(fMat_Z) / fMat_A) * (fDensity / 10.0) * kM_e; + + // Integrate dΩ²/dE = ω²_unit / |dEdx(E)|³ using trapezoidal rule on the same log-spaced grid. + std::vector rangeVar(nPoints, 0.0); + for (int i = 1; i < nPoints; ++i) { + double dE = energies[i] - energies[i - 1]; + double d0 = dedxValues[i - 1]; + double d1 = dedxValues[i]; + double f0 = (d0 > 0) ? 1.0 / (d0 * d0 * d0) : 0.0; + double f1 = (d1 > 0) ? 1.0 / (d1 * d1 * d1) : 0.0; + rangeVar[i] = rangeVar[i - 1] + omega2_unit * 0.5 * (f0 + f1) * dE; + } + fRangeVariance = tk::spline(energies, rangeVar); } double AtELossBetheBloch::GetdEdx_formula(double energy) const @@ -199,15 +217,12 @@ double AtELossBetheBloch::GetEnergy(double energyIni, double distance) const double AtELossBetheBloch::GetElossStraggling(double energyIni, double energyFin) const { - double dx_mm = GetRange(energyIni, energyFin); - if (dx_mm <= 0) + if (energyIni <= energyFin) return 0; - - // Bohr approximation: σ² = K·z²·(Z/A)·ρ·Δx[cm]·mₑc² - double z = IsElectron() ? 1.0 : fPart_q; - double dx_cm = dx_mm * 0.1; - double sigma2 = kK * z * z * (static_cast(fMat_Z) / fMat_A) * fDensity * dx_cm * kM_e; - return std::sqrt(sigma2); + double omega2 = GetRangeVariance(energyIni) - GetRangeVariance(energyFin); + if (omega2 <= 0) + return 0; + return std::abs(GetdEdx(energyIni)) * std::sqrt(omega2); } double AtELossBetheBloch::GetdEdxStraggling(double energyIni, double energyFin) const @@ -215,7 +230,19 @@ double AtELossBetheBloch::GetdEdxStraggling(double energyIni, double energyFin) double dx_mm = GetRange(energyIni, energyFin); if (dx_mm <= 0) return 0; - return GetElossStraggling(energyIni, energyFin) / dx_mm; + double omega2 = GetRangeVariance(energyIni) - GetRangeVariance(energyFin); + if (omega2 <= 0) + return 0; + return std::abs(GetdEdx(energyIni)) * std::sqrt(omega2) / dx_mm; +} + +double AtELossBetheBloch::GetRangeVariance(double energy) const +{ + if (energy <= fRangeVariance.get_x_min()) + return 0; + if (energy >= fRangeVariance.get_x_max()) + return fRangeVariance(fRangeVariance.get_x_max()); + return fRangeVariance(energy); } } // namespace AtTools diff --git a/AtTools/AtELossBetheBloch.h b/AtTools/AtELossBetheBloch.h index a2bcb515c..8dadd222e 100644 --- a/AtTools/AtELossBetheBloch.h +++ b/AtTools/AtELossBetheBloch.h @@ -28,7 +28,8 @@ class AtELossBetheBloch : public AtELossModel { int fMat_A; // target mass number (g/mol) double fI_MeV; // mean excitation energy in MeV - tk::spline fdXdE; // spline of dx/dE vs. energy; integral cached for O(log n) GetRange + tk::spline fdXdE; // spline of dx/dE vs. energy; integral cached for O(log n) GetRange + tk::spline fRangeVariance; // spline of Ω²(E) [mm²] vs. energy; accumulated Bohr range variance public: /** @@ -57,6 +58,7 @@ class AtELossBetheBloch : public AtELossModel { virtual double GetEnergy(double energyIni, double distance) const override; virtual double GetElossStraggling(double energyIni, double energyFin) const override; virtual double GetdEdxStraggling(double energyIni, double energyFin) const override; + virtual double GetRangeVariance(double energy) const override; private: bool IsElectron() const { return std::abs(fPart_mass - kM_e) < 0.01; } diff --git a/AtTools/AtELossBetheBlochTest.cxx b/AtTools/AtELossBetheBlochTest.cxx index 2dfcaec2e..d223e3bcc 100644 --- a/AtTools/AtELossBetheBlochTest.cxx +++ b/AtTools/AtELossBetheBlochTest.cxx @@ -145,15 +145,79 @@ TEST_F(AtELossBetheBlochFixture, SettersRebuildSpline) EXPECT_NE(m.GetdEdx(1.0), dedx_before); } -TEST_F(AtELossBetheBlochFixture, BohrStragglingSanity) +// Proton in Ar gas (same material as PiondEdx test). +class AtELossBetheBlochArFixture : public ::testing::Test { +protected: + AtELossBetheBloch model; + // Ar: Z=18, A=40, ρ=1.65e-3 g/cm³, I=188 eV + AtELossBetheBlochArFixture() : model(1.0, kProtonMass, 18, 40, 1.65e-3, 188.0) {} +}; + +TEST_F(AtELossBetheBlochArFixture, RangeVarianceZeroAtLowEnergy) +{ + EXPECT_DOUBLE_EQ(model.GetRangeVariance(0.0), 0.0); + EXPECT_DOUBLE_EQ(model.GetRangeVariance(1e-6), 0.0); +} + +TEST_F(AtELossBetheBlochArFixture, RangeVarianceMonotonicallyIncreasing) +{ + double prev = 0.0; + for (double E : {1.0, 2.0, 5.0, 10.0, 20.0, 50.0}) { + double rv = model.GetRangeVariance(E); + EXPECT_GT(rv, prev) << "RangeVariance not increasing at E=" << E; + prev = rv; + } +} + +TEST_F(AtELossBetheBlochArFixture, StragglingScalesWithSqrtRange) +{ + // For a factor-4 path length ratio, straggling should scale as sqrt(path) → factor ~2. + // Use two different path lengths by picking final energies via GetEnergy. + double E0 = 20.0; + double Ef1 = model.GetEnergy(E0, 500.0); + double Ef2 = model.GetEnergy(E0, 2000.0); // 4× longer path + ASSERT_GT(Ef1, 0.0); + ASSERT_GT(Ef2, 0.0); + + double sigma1 = model.GetElossStraggling(E0, Ef1); + double sigma2 = model.GetElossStraggling(E0, Ef2); + ASSERT_GT(sigma1, 0.0); + ASSERT_GT(sigma2, 0.0); + + // σ scales roughly as sqrt(Δx), so σ2/σ1 ≈ sqrt(4) = 2. + // Allow 20% tolerance: exact scaling only holds when dEdx is constant over the path. + double ratio = sigma2 / sigma1; + EXPECT_NEAR(ratio, 2.0, 0.20 * 2.0); +} + +TEST_F(AtELossBetheBlochArFixture, dEdxStragglingConsistency) { - double E0 = 5.0; - double Ef = model.GetEnergy(E0, 100.0); - EXPECT_GT(Ef, 0.0); + // GetdEdxStraggling * GetRange must equal GetElossStraggling (up to floating point). + double E0 = 10.0; + double Ef = model.GetEnergy(E0, 1000.0); + ASSERT_GT(Ef, 0.0); + + double sigmaE = model.GetElossStraggling(E0, Ef); + double dx = model.GetRange(E0, Ef); + double sigmaDedx = model.GetdEdxStraggling(E0, Ef); + + ASSERT_GT(sigmaE, 0.0); + ASSERT_GT(dx, 0.0); + EXPECT_NEAR(sigmaDedx * dx, sigmaE, 1e-9 * sigmaE); +} +TEST_F(AtELossBetheBlochArFixture, StragglingBroadRange) +{ + // Physical bound: straggling cannot exceed energy loss (σ(ΔE) < ΔE). + double E0 = 10.0; + double Ef = 5.0; + double dE = E0 - Ef; double sigma = model.GetElossStraggling(E0, Ef); + EXPECT_GT(sigma, 0.0); - EXPECT_LT(sigma, E0); // straggling must be less than total energy + EXPECT_LT(sigma, dE); + // Also check the ratio is sub-50% (physically meaningful bound) + EXPECT_LT(sigma / dE, 0.5); } TEST_F(AtELossBetheBlochFixture, BlochApprox) From e0b2ba242b1d6f8aa05d01780e7de66737b9d4ea Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Tue, 31 Mar 2026 20:18:02 -0400 Subject: [PATCH 02/10] Fix Bethe-Bloch straggling and add CATIMA benchmark tests MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit - GetElossStraggling/GetdEdxStraggling: evaluate dEdx at energyFin (Tschalar formula requires stopping power at the exit point, not entrance; error was largest near the Bragg peak) - BuildSpline: replace manual trapezoidal rule for Ω²(E) with a spline of 1/dEdx³ and segment-wise Simpson integration (exact for cubics) - AtELossCATIMA: add fConfig member and SetConfig() so callers can tune the CATIMA calculation model; thread fConfig through all catima API calls (default_config preserves existing behaviour) - AtELossBetheBlochTest: add AtELossBetheBlochVsCATIMAFixture comparing GetRangeVariance and GetElossStraggling against pure-Bohr CATIMA (z_eff_type::none) at narrow and wide energy traversals up to 10→1 MeV --- AtTools/AtELossBetheBloch.cxx | 22 +++++----- AtTools/AtELossBetheBlochTest.cxx | 68 +++++++++++++++++++++++++++++++ AtTools/AtELossCATIMA.cxx | 14 +++---- AtTools/AtELossCATIMA.h | 3 ++ 4 files changed, 89 insertions(+), 18 deletions(-) diff --git a/AtTools/AtELossBetheBloch.cxx b/AtTools/AtELossBetheBloch.cxx index 90437c094..4d36362d4 100644 --- a/AtTools/AtELossBetheBloch.cxx +++ b/AtTools/AtELossBetheBloch.cxx @@ -85,16 +85,16 @@ void AtELossBetheBloch::BuildSpline(double E_min_MeV, double E_max_MeV, int nPoi double z = IsElectron() ? 1.0 : fPart_q; double omega2_unit = kK * z * z * (static_cast(fMat_Z) / fMat_A) * (fDensity / 10.0) * kM_e; - // Integrate dΩ²/dE = ω²_unit / |dEdx(E)|³ using trapezoidal rule on the same log-spaced grid. + // Build spline of the integrand dΩ²/dE = 1/|dEdx|³, then integrate segment-by-segment using + // Simpson's rule (exact for cubics, O(h⁴)) rather than the O(h²) trapezoidal rule. + std::vector integrandValues(nPoints); + for (int i = 0; i < nPoints; ++i) + integrandValues[i] = (dedxValues[i] > 0) ? 1.0 / (dedxValues[i] * dedxValues[i] * dedxValues[i]) : 0.0; + tk::spline integrand(energies, integrandValues); + std::vector rangeVar(nPoints, 0.0); - for (int i = 1; i < nPoints; ++i) { - double dE = energies[i] - energies[i - 1]; - double d0 = dedxValues[i - 1]; - double d1 = dedxValues[i]; - double f0 = (d0 > 0) ? 1.0 / (d0 * d0 * d0) : 0.0; - double f1 = (d1 > 0) ? 1.0 / (d1 * d1 * d1) : 0.0; - rangeVar[i] = rangeVar[i - 1] + omega2_unit * 0.5 * (f0 + f1) * dE; - } + for (int i = 1; i < nPoints; ++i) + rangeVar[i] = rangeVar[i - 1] + omega2_unit * integrand.integrate(energies[i - 1], energies[i]); fRangeVariance = tk::spline(energies, rangeVar); } @@ -222,7 +222,7 @@ double AtELossBetheBloch::GetElossStraggling(double energyIni, double energyFin) double omega2 = GetRangeVariance(energyIni) - GetRangeVariance(energyFin); if (omega2 <= 0) return 0; - return std::abs(GetdEdx(energyIni)) * std::sqrt(omega2); + return std::abs(GetdEdx(energyFin)) * std::sqrt(omega2); } double AtELossBetheBloch::GetdEdxStraggling(double energyIni, double energyFin) const @@ -233,7 +233,7 @@ double AtELossBetheBloch::GetdEdxStraggling(double energyIni, double energyFin) double omega2 = GetRangeVariance(energyIni) - GetRangeVariance(energyFin); if (omega2 <= 0) return 0; - return std::abs(GetdEdx(energyIni)) * std::sqrt(omega2) / dx_mm; + return std::abs(GetdEdx(energyFin)) * std::sqrt(omega2) / dx_mm; } double AtELossBetheBloch::GetRangeVariance(double energy) const diff --git a/AtTools/AtELossBetheBlochTest.cxx b/AtTools/AtELossBetheBlochTest.cxx index d223e3bcc..f1979fa75 100644 --- a/AtTools/AtELossBetheBlochTest.cxx +++ b/AtTools/AtELossBetheBlochTest.cxx @@ -1,5 +1,8 @@ #include "AtELossBetheBloch.h" +#include "AtELossCATIMA.h" +#include +#include #include #include @@ -232,3 +235,68 @@ TEST_F(AtELossBetheBlochFixture, BlochApprox) // Bloch approximation should agree to within 30% of the explicit I value EXPECT_NEAR(rangeBloch, rangeExact, 0.30 * rangeExact); } + +/** + * Benchmark fixture: compare AtELossBetheBloch against AtELossCATIMA configured for pure Bohr + * straggling (z_eff_type::none = bare charge, default calculation = bohr mode). + * + * System: proton in H₂ at 600 Torr (ρ = 6.5643e-5 g/cm³), identical to AtELossCATIMATestFixture. + * Expected agreement: within 10% (residual difference from Lindhard X × γ² ≈ 1 for protons at 1–10 MeV). + */ +class AtELossBetheBlochVsCATIMAFixture : public ::testing::Test { +protected: + static constexpr double kProtonMassAmu = 1.007825031898; + static constexpr double kDensity = 6.5643e-5; // g/cm³ + + AtTools::AtELossBetheBloch bb; + AtTools::AtELossCATIMA catima; + + AtELossBetheBlochVsCATIMAFixture() + : bb(1.0, kProtonMass, 1, 1, kDensity, 19.2), catima(kDensity, catima::Material(1, 1)) + { + catima.SetProjectile(1, 1, kProtonMassAmu); + + // Configure CATIMA to pure Bohr: bare charge (no effective-charge model), + // default calculation=bohr already disables the ATIMA correction. + catima::Config pureBohrCfg; + pureBohrCfg.z_effective = catima::z_eff_type::none; + catima.SetConfig(pureBohrCfg); + } +}; + +TEST_F(AtELossBetheBlochVsCATIMAFixture, CATIMAComparison_RangeVariance) +{ + // Compare accumulated Bohr range variance Ω²(E) [mm²] at three energies. + // Both models implement the same Bohr formula; agreement within 10% is expected. + for (double E : {1.0, 5.0, 10.0}) { + double bbVal = bb.GetRangeVariance(E); + double catimaVal = catima.GetRangeVariance(E); + ASSERT_GT(bbVal, 0.0) << "BB range variance zero at E=" << E; + ASSERT_GT(catimaVal, 0.0) << "CATIMA range variance zero at E=" << E; + EXPECT_NEAR(bbVal, catimaVal, 0.10 * catimaVal) << "Range variance mismatch at E=" << E << " MeV"; + } +} + +TEST_F(AtELossBetheBlochVsCATIMAFixture, CATIMAComparison_EnergyStraggling) +{ + // Reference energy pairs from AtELossCATIMATestFixture::TestEnergyLossStraggling. + // CATIMA (pure Bohr config) is used as the reference; BB must agree within 10%. + const double mass = kProtonMassAmu; + struct Case { + double eIni, eFin; + }; + const Case cases[] = { + {1.0, 0.75 * mass}, // narrow: ~25% energy loss + {5.0, 3.58164 * mass}, // narrow: ~28% energy loss + {5.0, 1.0}, // wide: 80% energy loss + {10.0, 1.0}, // wide: 90% energy loss, exercises Bragg-peak region + }; + for (auto &c : cases) { + double bbSigma = bb.GetElossStraggling(c.eIni, c.eFin); + double catimaSigma = catima.GetElossStraggling(c.eIni, c.eFin); + ASSERT_GT(bbSigma, 0.0) << "BB straggling zero for E0=" << c.eIni; + ASSERT_GT(catimaSigma, 0.0) << "CATIMA straggling zero for E0=" << c.eIni; + EXPECT_NEAR(bbSigma, catimaSigma, 0.10 * catimaSigma) + << "Energy straggling mismatch for E0=" << c.eIni << " MeV"; + } +} diff --git a/AtTools/AtELossCATIMA.cxx b/AtTools/AtELossCATIMA.cxx index 3f85a981d..e8b7919d2 100644 --- a/AtTools/AtELossCATIMA.cxx +++ b/AtTools/AtELossCATIMA.cxx @@ -23,7 +23,7 @@ double AtELossCATIMA::GetdEdx(double energy) const return 0; } - catima::Result result = catima::calculate(*fProjectile, *fMaterial, energy / fProjectileMassAmu); + catima::Result result = catima::calculate(*fProjectile, *fMaterial, energy / fProjectileMassAmu, fConfig); double dEdx = result.dEdxi * fDensity; // MeV/cm return dEdx / 10.0; // convert to MeV/mm } @@ -38,13 +38,13 @@ double AtELossCATIMA::GetRange(double energyIni, double energyFin) const if (energyFin == 0) { fProjectile->T = energyIni / fProjectileMassAmu; - return catima::range(*fProjectile, *fMaterial) / fDensity * 10.; + return catima::range(*fProjectile, *fMaterial, fConfig) / fDensity * 10.; } double remainingEnergy{energyIni}; double range{0}; while (remainingEnergy > energyFin) { - catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu); + catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); double dEdx = result.dEdxi * fDensity; double DE = dEdx * fRangeStepSize / 10.; @@ -65,7 +65,7 @@ double AtELossCATIMA::GetEnergy(double energyIni, double distance) const double remainingEnergy{energyIni}; double range{0}; while (range < distance) { - catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu); + catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); double dEdx = result.dEdxi * fDensity; double DE{}; @@ -90,7 +90,7 @@ double AtELossCATIMA::GetRangeVariance(double energy) const return 0; } auto range_var = - catima::range_variance(*fProjectile, energy / fProjectileMassAmu, *fMaterial); // range var in (g/cm^2)^2 + catima::range_variance(*fProjectile, energy / fProjectileMassAmu, *fMaterial, fConfig); // range var in (g/cm^2)^2 LOG(debug) << "Range variance in (g/cm^2)^2: " << range_var << " for energy: " << energy / fProjectileMassAmu << " MeV/u"; range_var /= fDensity * fDensity; // convert to (cm)^2 @@ -110,7 +110,7 @@ double AtELossCATIMA::GetElossStraggling(double energyIni, double energyFin) con return 0; } auto energy_strag = catima::energy_straggling_from_E(*fProjectile, energyIni / fProjectileMassAmu, - energyFin / fProjectileMassAmu, *fMaterial); + energyFin / fProjectileMassAmu, *fMaterial, fConfig); return energy_strag; } double AtELossCATIMA::GetdEdxStraggling(double energyIni, double energyFin) const @@ -142,7 +142,7 @@ AtELossCATIMA::GetBraggCurve(double energy, double rangeStepSize, double totalFr double range{}; while (remainingEnergy / energy > totalFractionELoss) { - catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu); + catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); double dEdx = result.dEdxi * fDensity; braggCurve.push_back(std::make_pair(dEdx, range)); diff --git a/AtTools/AtELossCATIMA.h b/AtTools/AtELossCATIMA.h index bf63862e7..f7b3fbb63 100644 --- a/AtTools/AtELossCATIMA.h +++ b/AtTools/AtELossCATIMA.h @@ -5,6 +5,7 @@ #include "AtELossModel.h" #include +#include #include #include #include @@ -19,6 +20,7 @@ class AtELossCATIMA : public AtELossModel { double fProjectileMassAmu{-1}; /// Mass of the projectile in amu (atomic mass units). double fRangeStepSize{0.1}; // mm + catima::Config fConfig{catima::default_config}; public: /** @@ -68,6 +70,7 @@ class AtELossCATIMA : public AtELossModel { * @param[in] stepSize The step size used for ranges. It must be input in mm. */ void SetRangeStepSize(double stepSize) { fRangeStepSize = stepSize; } + void SetConfig(catima::Config cfg) { fConfig = cfg; } }; } // namespace AtTools From ae431bbf1e5494b8c1846907f60e3ab669111d44 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Tue, 31 Mar 2026 20:45:52 -0400 Subject: [PATCH 03/10] Add energy loss model comparison macro Adds macro/tests/ELossComparison/ELossComparison.C with four scenarios (proton/deuteron/alpha in D2, proton in He) comparing BetheBloch, CATIMA, and SRIM Table models via dE/dx, range, Bragg curve, and straggling plots. --- .gitignore | 1 + macro/tests/ELossComparison/ELossComparison.C | 418 ++++++++++++++++++ macro/tests/ELossComparison/runTests.sh | 18 + 3 files changed, 437 insertions(+) create mode 100644 macro/tests/ELossComparison/ELossComparison.C create mode 100755 macro/tests/ELossComparison/runTests.sh diff --git a/.gitignore b/.gitignore index f0c9c9839..709e805eb 100755 --- a/.gitignore +++ b/.gitignore @@ -1,3 +1,4 @@ +*.pdf *.swp *.root *.log diff --git a/macro/tests/ELossComparison/ELossComparison.C b/macro/tests/ELossComparison/ELossComparison.C new file mode 100644 index 000000000..ddc91a742 --- /dev/null +++ b/macro/tests/ELossComparison/ELossComparison.C @@ -0,0 +1,418 @@ +/** + * ELossComparison.C + * + * Comparative plots of the three energy loss models in ATTPCROOTv2: + * - AtELossBetheBloch (analytic PDG Bethe-Bloch formula) + * - AtELossCATIMA (CATIMA library wrapper) + * - AtELossTable (SRIM lookup table with spline interpolation) + * + * Four scenarios are plotted: + * 1. Proton in H2 gas @ 600 Torr (all three models) + * 2. Deuteron in D2 gas @ 600 Torr (BB + CATIMA + Table) + * 3. Alpha (4He) in D2 gas @ 600 Torr (BB + CATIMA + Table) + * 4. Proton in He gas @ 700 Torr (BB + CATIMA + Table) + * + * Each scenario produces a canvas with four pads: + * Top-left dE/dx vs kinetic energy (log-log) + * Top-right Range vs kinetic energy (log-log) + * Bottom-left Bragg curve (dE/dx vs distance from entry) + * Bottom-right Range straggling sigma vs kinetic energy (log-log) + * Note: Table straggling requires a full SRIM file with the + * "Multiply" conversion section. Simple SRIM tables (no + * conversion header) are loaded without straggling data and + * are omitted from the straggling pad. + * + * Usage (after sourcing build/config.sh): + * root -l -b -q ELossComparison.C + * + * Output: ELossComparison.pdf (four pages, one per scenario) + */ + +// --------------------------------------------------------------------------- +// Unit conversion helpers (energy -> MeV, length -> mm) +// --------------------------------------------------------------------------- +double UnitToMeV(const std::string &u) +{ + if (u == "eV") + return 1e-6; + if (u == "keV") + return 1e-3; + if (u == "MeV") + return 1.0; + if (u == "GeV") + return 1e3; + return 0.0; +} + +double UnitToMM(const std::string &u) +{ + if (u == "um") + return 1e-3; + if (u == "mm") + return 1.0; + if (u == "cm") + return 10.0; + if (u == "m") + return 1000.0; + return 0.0; +} + +// --------------------------------------------------------------------------- +// LoadSimpleSrimTable +// +// Load a SRIM-format stopping table that has the standard column layout but +// lacks the "Target Density" header and "Multiply" conversion section. +// +// SRIM stores dE/dx in MeV/(mg/cm2). With known density rho (g/cm3): +// dEdx [MeV/mm] = dEdx [MeV/(mg/cm2)] x rho [g/cm3] x 100 +// +// Range straggling cannot be loaded because AtELossTable::LoadRangeVariance +// is private; the returned model has no straggling data. +// --------------------------------------------------------------------------- +AtTools::AtELossTable *LoadSimpleSrimTable(const std::string &fileName, double density) +{ + std::ifstream file(fileName); + if (!file.is_open()) { + ::Error("LoadSimpleSrimTable", "Cannot open file: %s", fileName.c_str()); + return nullptr; + } + + const double conversion = density * 100.0; // MeV/(mg/cm2) -> MeV/mm + + std::vector energy; + std::vector dEdX; + + std::string line; + while (std::getline(file, line)) { + std::istringstream iss(line); + std::vector tok; + std::string t; + while (iss >> t) + tok.push_back(t); + + // Need: energy Eunit dEdx_e dEdx_n range Runit strag Sunit (>=8 tokens) + if (tok.size() < 8) + continue; + + double en, dedxE, dedxN; + try { + en = std::stod(tok[0]); + dedxE = std::stod(tok[2]); + dedxN = std::stod(tok[3]); + } catch (...) { + continue; + } + + double enMeV = en * UnitToMeV(tok[1]); + double dedxMeVmm = (dedxE + dedxN) * conversion; + + if (enMeV <= 0.0 || dedxMeVmm <= 0.0) + continue; + + energy.push_back(enMeV); + dEdX.push_back(dedxMeVmm); + } + + if (energy.empty()) { + ::Error("LoadSimpleSrimTable", "No valid rows parsed from %s", fileName.c_str()); + return nullptr; + } + + return new AtTools::AtELossTable(energy, dEdX, density); +} + +// --------------------------------------------------------------------------- +// Graph-filling helpers: sample a model over a log-spaced energy grid +// --------------------------------------------------------------------------- +TGraph *MakeDedxGraph(AtTools::AtELossModel *m, double eMin, double eMax, int nPts = 200) +{ + auto *g = new TGraph(nPts); + for (int i = 0; i < nPts; ++i) { + double e = eMin * TMath::Power(eMax / eMin, double(i) / (nPts - 1)); + g->SetPoint(i, e, m->GetdEdx(e)); + } + return g; +} + +TGraph *MakeRangeGraph(AtTools::AtELossModel *m, double eMin, double eMax, int nPts = 200) +{ + auto *g = new TGraph(nPts); + for (int i = 0; i < nPts; ++i) { + double e = eMin * TMath::Power(eMax / eMin, double(i) / (nPts - 1)); + g->SetPoint(i, e, m->GetRange(e)); + } + return g; +} + +// hasStrag: pass false for table models loaded without straggling data +TGraph *MakeStragGraph(AtTools::AtELossModel *m, double eMin, double eMax, bool hasStrag, int nPts = 200) +{ + auto *g = new TGraph(); + if (!hasStrag) + return g; + for (int i = 0; i < nPts; ++i) { + double e = eMin * TMath::Power(eMax / eMin, double(i) / (nPts - 1)); + double s = m->GetRangeStraggling(e); + if (s > 0.0) + g->AddPoint(e, s); + } + return g; +} + +TGraph *MakeBraggGraph(AtTools::AtELossModel *m, double energy, double stepMM = 0.5) +{ + auto pts = m->GetBraggCurve(energy, stepMM); + auto *g = new TGraph(int(pts.size())); + for (int i = 0; i < int(pts.size()); ++i) + g->SetPoint(i, pts[i].second, pts[i].first); + return g; +} + +// --------------------------------------------------------------------------- +// ModelEntry: bundles a model pointer with display properties +// --------------------------------------------------------------------------- +struct ModelEntry { + AtTools::AtELossModel *model; + std::string name; + int color; + int lineStyle; + bool hasStraggling; +}; + +// --------------------------------------------------------------------------- +// DrawScenario: fills a pre-divided 4-pad canvas for one particle/material +// --------------------------------------------------------------------------- +void DrawScenario(const std::string &title, const std::vector &entries, double eMin, double eMax, + double braggEnergy, TCanvas *c) +{ + c->SetTitle(title.c_str()); + c->Divide(2, 2, 0.005, 0.005); + + auto StyleGraph = [](TGraph *g, int col, int sty) { + g->SetLineColor(col); + g->SetLineStyle(sty); + g->SetLineWidth(2); + }; + + // --- Pad 1: dE/dx vs E --- + { + c->cd(1)->SetLogx(); + c->cd(1)->SetLogy(); + c->cd(1)->SetGrid(); + auto *mg = new TMultiGraph(); + auto *leg = new TLegend(0.55, 0.65, 0.88, 0.88); + leg->SetBorderSize(0); + for (auto &e : entries) { + auto *g = MakeDedxGraph(e.model, eMin, eMax); + StyleGraph(g, e.color, e.lineStyle); + mg->Add(g, "L"); + leg->AddEntry(g, e.name.c_str(), "L"); + } + mg->SetTitle(Form("%s;Kinetic energy (MeV);dE/dx (MeV/mm)", title.c_str())); + mg->Draw("A"); + leg->Draw(); + } + + // --- Pad 2: Range vs E --- + { + c->cd(2)->SetLogx(); + c->cd(2)->SetLogy(); + c->cd(2)->SetGrid(); + auto *mg = new TMultiGraph(); + auto *leg = new TLegend(0.15, 0.65, 0.48, 0.88); + leg->SetBorderSize(0); + for (auto &e : entries) { + auto *g = MakeRangeGraph(e.model, eMin, eMax); + StyleGraph(g, e.color, e.lineStyle); + mg->Add(g, "L"); + leg->AddEntry(g, e.name.c_str(), "L"); + } + mg->SetTitle(";Kinetic energy (MeV);Range (mm)"); + mg->Draw("A"); + leg->Draw(); + } + + // --- Pad 3: Bragg curve --- + { + c->cd(3)->SetGrid(); + auto *mg = new TMultiGraph(); + auto *leg = new TLegend(0.15, 0.65, 0.48, 0.88); + leg->SetBorderSize(0); + for (auto &e : entries) { + auto *g = MakeBraggGraph(e.model, braggEnergy); + StyleGraph(g, e.color, e.lineStyle); + mg->Add(g, "L"); + leg->AddEntry(g, e.name.c_str(), "L"); + } + mg->SetTitle(Form("Bragg curve (E_{0} = %.1f MeV);Distance from entry (mm);dE/dx (MeV/mm)", braggEnergy)); + mg->Draw("A"); + leg->Draw(); + } + + // --- Pad 4: Range straggling vs E --- + { + c->cd(4)->SetLogx(); + c->cd(4)->SetLogy(); + c->cd(4)->SetGrid(); + auto *mg = new TMultiGraph(); + auto *leg = new TLegend(0.15, 0.65, 0.48, 0.88); + leg->SetBorderSize(0); + bool anyDrawn = false; + for (auto &e : entries) { + auto *g = MakeStragGraph(e.model, eMin, eMax, e.hasStraggling); + if (g->GetN() == 0) { + delete g; + continue; + } + StyleGraph(g, e.color, e.lineStyle); + mg->Add(g, "L"); + leg->AddEntry(g, e.name.c_str(), "L"); + anyDrawn = true; + } + if (anyDrawn) { + mg->SetTitle(";Kinetic energy (MeV);Range straggling #sigma (mm)"); + mg->Draw("A"); + } + leg->Draw(); + } +} + +// =========================================================================== +// Main macro +// =========================================================================== +void ELossComparison() +{ + TString vmcDir = getenv("VMCWORKDIR") ? getenv("VMCWORKDIR") : "../../.."; + TString resDir = vmcDir + "/resources/energy_loss"; + + // ------------------------------------------------------------------------- + // Physical parameters + // ------------------------------------------------------------------------- + + // Particle rest masses (MeV/c2) for Bethe-Bloch + const double kMproton = 938.272; + const double kMdeuteron = 1875.61; + const double kMalpha = 3727.38; + + // Particle masses in amu for CATIMA + const double kMproton_amu = 1.00728; + const double kMdeuteron_amu = 2.01355; + const double kMalpha_amu = 4.00151; + + // Gas densities (g/cm3) at stated pressures, room temperature (~20 C) + const double kRho_H2_600 = 6.5643e-5; // H2 @ 600 Torr (CATIMA/LISE reference tests) + const double kRho_D2_600 = 1.313e-4; // D2 @ 600 Torr (scaled by molar-mass ratio from H2) + const double kRho_He_700 = 1.53e-4; // He @ 700 Torr (He STP x 700/760 x 273/293) + + // Mean excitation energies (eV) -- PDG 2022 Table 34.1 + const double kI_H = 19.2; + const double kI_He = 41.8; + + // ========================================================================= + // Scenario 1: Proton in H2 @ 600 Torr + // All three models. HinH.txt has the full SRIM format (density + Multiply) + // so the table model includes straggling data. + // ========================================================================= + auto bb_p_H2 = new AtTools::AtELossBetheBloch(1.0, kMproton, 1, 1, kRho_H2_600, kI_H); + + auto catima_p_H2 = + new AtTools::AtELossCATIMA(kRho_H2_600, std::vector>{{1, 1, 1}}); + catima_p_H2->SetProjectile(1, 1, kMproton_amu); + + auto table_p_H2 = new AtTools::AtELossTable(0); + table_p_H2->LoadSrimTable(Form("%s/HinH.txt", resDir.Data())); + table_p_H2->SetDensity(kRho_H2_600); + + // ========================================================================= + // Scenario 2: Deuteron in D2 @ 600 Torr + // ========================================================================= + auto bb_d_D2 = new AtTools::AtELossBetheBloch(1.0, kMdeuteron, 1, 2, kRho_D2_600, kI_H); + + auto catima_d_D2 = + new AtTools::AtELossCATIMA(kRho_D2_600, std::vector>{{2, 1, 1}}); + catima_d_D2->SetProjectile(2, 1, kMdeuteron_amu); + + auto table_d_D2 = LoadSimpleSrimTable(Form("%s/deuteron_D2_600torr.txt", resDir.Data()), kRho_D2_600); + + // ========================================================================= + // Scenario 3: Alpha (4He) in D2 @ 600 Torr + // ========================================================================= + auto bb_a_D2 = new AtTools::AtELossBetheBloch(2.0, kMalpha, 1, 2, kRho_D2_600, kI_H); + + auto catima_a_D2 = + new AtTools::AtELossCATIMA(kRho_D2_600, std::vector>{{2, 1, 1}}); + catima_a_D2->SetProjectile(4, 2, kMalpha_amu); + + auto table_a_D2 = LoadSimpleSrimTable(Form("%s/alpha_D2_600torr.txt", resDir.Data()), kRho_D2_600); + + // ========================================================================= + // Scenario 4: Proton in He @ 700 Torr + // ========================================================================= + auto bb_p_He = new AtTools::AtELossBetheBloch(1.0, kMproton, 2, 4, kRho_He_700, kI_He); + + auto catima_p_He = + new AtTools::AtELossCATIMA(kRho_He_700, std::vector>{{4, 2, 1}}); + catima_p_He->SetProjectile(1, 1, kMproton_amu); + + auto table_p_He = LoadSimpleSrimTable(Form("%s/proton_He_700torr.txt", resDir.Data()), kRho_He_700); + + // ========================================================================= + // Color / line-style scheme + // ========================================================================= + const int kColBB = kRed + 1; + const int kColCATIMA = kBlue + 1; + const int kColTable = kGreen + 2; + + // ========================================================================= + // Draw and save + // ========================================================================= + auto *c1 = new TCanvas("c1", "Proton in H2", 1200, 900); + { + std::vector ents = { + {bb_p_H2, "Bethe-Bloch", kColBB, 1, true}, + {catima_p_H2, "CATIMA", kColCATIMA, 2, true}, + {table_p_H2, "SRIM Table", kColTable, 3, true}, + }; + DrawScenario("Proton in H_{2} (600 Torr)", ents, 0.1, 10.0, 5.0, c1); + } + + auto *c2 = new TCanvas("c2", "Deuteron in D2", 1200, 900); + { + std::vector ents = { + {bb_d_D2, "Bethe-Bloch", kColBB, 1, true}, + {catima_d_D2, "CATIMA", kColCATIMA, 2, true}, + }; + if (table_d_D2) + ents.push_back({table_d_D2, "SRIM Table", kColTable, 3, false}); + DrawScenario("Deuteron in D_{2} (600 Torr)", ents, 0.2, 20.0, 10.0, c2); + } + + auto *c3 = new TCanvas("c3", "Alpha in D2", 1200, 900); + { + std::vector ents = { + {bb_a_D2, "Bethe-Bloch", kColBB, 1, true}, + {catima_a_D2, "CATIMA", kColCATIMA, 2, true}, + }; + if (table_a_D2) + ents.push_back({table_a_D2, "SRIM Table", kColTable, 3, false}); + DrawScenario("Alpha (^{4}He) in D_{2} (600 Torr)", ents, 0.5, 30.0, 10.0, c3); + } + + auto *c4 = new TCanvas("c4", "Proton in He", 1200, 900); + { + std::vector ents = { + {bb_p_He, "Bethe-Bloch", kColBB, 1, true}, + {catima_p_He, "CATIMA", kColCATIMA, 2, true}, + }; + if (table_p_He) + ents.push_back({table_p_He, "SRIM Table", kColTable, 3, false}); + DrawScenario("Proton in ^{4}He (700 Torr)", ents, 0.1, 10.0, 5.0, c4); + } + + c1->Print("ELossComparison.pdf("); + c2->Print("ELossComparison.pdf"); + c3->Print("ELossComparison.pdf"); + c4->Print("ELossComparison.pdf)"); + + ::Info("ELossComparison", "Saved to ELossComparison.pdf"); +} diff --git a/macro/tests/ELossComparison/runTests.sh b/macro/tests/ELossComparison/runTests.sh new file mode 100755 index 000000000..04bf95d22 --- /dev/null +++ b/macro/tests/ELossComparison/runTests.sh @@ -0,0 +1,18 @@ +#!/bin/bash +# Run the ELossComparison integration macro. +# Must be called after sourcing build/config.sh (sets VMCWORKDIR, LD_LIBRARY_PATH, etc.) + +set -e + +SCRIPT_DIR="$( cd "$( dirname "${BASH_SOURCE[0]}" )" && pwd )" +cd "${SCRIPT_DIR}" + +if [ -z "${VMCWORKDIR}" ]; then + echo "ERROR: VMCWORKDIR is not set. Source build/config.sh first." >&2 + exit 1 +fi + +echo "Running ELossComparison..." +root -l -b -q ELossComparison.C + +echo "Done. Output: ${SCRIPT_DIR}/ELossComparison.pdf" From 091dca258751bbc01b72ec41e033691618036814 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Tue, 31 Mar 2026 21:25:59 -0400 Subject: [PATCH 04/10] Add documentation useful for agents/humans Add documentation into four subdirectories directly under docs/: - tooling/: installation, daily-use (merged setup+building), testing - contributing/: contributor guide, new-module, code-style - reference/: modules, data-model, branch-io-contracts, macro-cookbook - subsystems/: simulation/reconstruction pipelines, generators, psa, energy-loss docs/development/ is reserved for branch-specific in-progress notes (fitting-status.md, UKF.md, etc.). Update all cross-links and the CLAUDE.md to point to the new docs rather than repeat the information. --- ...is-implementing-a-kalman-filter-for-us.txt | 867 ++++++++++++++++++ .claude/CLAUDE.md | 170 +--- .claude/plans/keen-marinating-wind.md | 295 ++++++ .claude/resume-ukf-task-integration.md | 38 + AGENTS.md | 57 ++ docs/contributing/code-style.md | 35 + docs/contributing/guide.md | 94 ++ docs/contributing/new-module.md | 71 ++ docs/development/fitting-status.md | 19 + docs/index.md | 34 + docs/reference/branch-io-contracts.md | 33 + docs/reference/data-model.md | 46 + docs/reference/macro-cookbook.md | 26 + docs/reference/modules.md | 28 + docs/subsystems/energy-loss.md | 59 ++ docs/subsystems/generators.md | 36 + docs/subsystems/psa.md | 40 + docs/subsystems/reconstruction-pipeline.md | 62 ++ docs/subsystems/simulation-pipeline.md | 54 ++ docs/tooling/daily-use.md | 61 ++ docs/tooling/installation.md | 50 + docs/tooling/testing.md | 59 ++ 22 files changed, 2106 insertions(+), 128 deletions(-) create mode 100644 .claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt create mode 100644 .claude/plans/keen-marinating-wind.md create mode 100644 .claude/resume-ukf-task-integration.md create mode 100644 AGENTS.md create mode 100644 docs/contributing/code-style.md create mode 100644 docs/contributing/guide.md create mode 100644 docs/contributing/new-module.md create mode 100644 docs/development/fitting-status.md create mode 100644 docs/index.md create mode 100644 docs/reference/branch-io-contracts.md create mode 100644 docs/reference/data-model.md create mode 100644 docs/reference/macro-cookbook.md create mode 100644 docs/reference/modules.md create mode 100644 docs/subsystems/energy-loss.md create mode 100644 docs/subsystems/generators.md create mode 100644 docs/subsystems/psa.md create mode 100644 docs/subsystems/reconstruction-pipeline.md create mode 100644 docs/subsystems/simulation-pipeline.md create mode 100644 docs/tooling/daily-use.md create mode 100644 docs/tooling/installation.md create mode 100644 docs/tooling/testing.md diff --git a/.claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt b/.claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt new file mode 100644 index 000000000..cf1499347 --- /dev/null +++ b/.claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt @@ -0,0 +1,867 @@ + +╭─── Claude Code v2.1.63 ──────────────────────────────────────────────────────╮ +│ │ Tips for getting │ +│ Welcome back Adam! │ started │ +│ │ Run /init to create a … │ +│ ▐▛███▜▌ │ ─────────────────────── │ +│ ▝▜█████▛▘ │ Recent activity │ +│ ▘▘ ▝▝ │ No recent activity │ +│ │ │ +│ Sonnet 4.6 · Claude Pro · aanthony3@gmail.com's │ │ +│ Organization │ │ +│ ~/fair_install/ATTPCROOTv2 │ │ +╰──────────────────────────────────────────────────────────────────────────────╯ + +❯ This branch is implementing a Kalman Filter for use in TPCs. You can see + tests showing how it is being used to fit tracks in @macro/tests/UKF/. + The class itself is located in + @AtReconstruction/AtFitter/OpenKF/kalman_filter/TrackFitterUKF* and it + depends on the propagator defined in @AtTools/AtPropagator*. This fitter + needs to be integrated into the task system of the codebase. The output + will be an AtFittedTrack which has a proposed class in the PR #254 which + is to merge the branch RealAurio:fitterRefactor. Please make a plan to + finalize this refactor and implement a task that will run the kalman + filter on tracks. The AtFitter task should share a base with the other + fitter tasks. At least the GenFit version, but also the MCFitter if it + combines well. + ⎿  Listed directory macro/tests/UKF/ + +● 3 Explore agents finished (ctrl+o to expand) + ├─ Explore UKF fitter and propagator classes · 12 tool uses · 55.6k tokens + │ ⎿  Done + ├─ Explore existing fitter tasks and task system · 25 tool uses · 50.2k + tokens + │ ⎿  Done + └─ Explore AtFittedTrack and PR #254 branch · 43 tool uses · 55.1k tokens + ⎿  Done + +● Let me read the critical files directly before designing the plan. + +● Read 6 files (ctrl+o to expand) + +● Now let me look at the CMakeLists files and the existing OpenKFTest to + understand build system integration: + +● Searched for 2 patterns, read 1 file (ctrl+o to expand) + +● Let me check the AtHitCluster and AtTrack data to understand what information + is available for UKF: + +● Searched for 2 patterns, read 1 file (ctrl+o to expand) + +● Now let me look at TrackFitterUKF's physics-specific interface in detail: + +● Searched for 1 pattern, read 1 file (ctrl+o to expand) + +● Now let me quickly check the ELoss model setup and check what the AtGenfit + does for momentum estimation: + +● Searched for 4 patterns (ctrl+o to expand) + +● Before finalizing the plan, I need to clarify a few key decisions: + ⎿  Tool use rejected with user message: The user wants to clarify these + questions. + This means they may have additional information, context or questions + for you. + Take their response into account and then reformulate the questions if + appropriate. + Start by asking them what they would like to clarify. + + Questions asked: + - "The fitterRefactor PR (#254) changes ~100+ files (namespace rename + AtFITTER→EventFit, AtFitter interface change from ProcessTracks() to + FitEvent()/GetFittedTrack(), new AtFittedTrack nested structs, removes PRA + data). How much of this refactor should we adopt?" + (No answer provided) + - "The UKF needs an initial momentum magnitude as a seed. What is the + primary strategy for providing this?" + (No answer provided) + - "What UKF output should be stored in AtFittedTrack beyond the standard + vertex state (energy, theta, phi, position)?" + (No answer provided) + +● What would you like to clarify about the questions I asked? + +✻ Baked for 4m 23s + +❯ I think you are comparing the wrong things. main is not the main branch, + develop is. Main is very old and should be ignored. That PR only touches + 36 files, not 100. Those changes are all important to understand to + construct a plan against them. To set the initial momentum, we will use a + circle fit like we do in the getfit case (PRA). Give me your options for + the last question again. + +● Good point - let me look at the fitterRefactor changes properly (against + develop) and how GenFit uses the circle fit for initial momentum. + +● Bash(git log --oneline RealAurio/fitterRefactor ^develop | head -20) + ⎿  11b7a6f9 Removed stuff from AtFitter that ended up being too case spec + ific to make sense to keep in the base class. Also fixed some serious + bugs that I found while actually implementing a fitter class in the br + … +16 lines (ctrl+o to expand) + +● Bash(git diff --name-only RealAurio/fitterRefactor develop -- 2>/dev/null | + head -50) + ⎿  AtData/AtDataLinkDef.h + AtData/AtFitMetadata.cxx + AtData/AtFitMetadata.h + … +27 lines (ctrl+o to expand) + +● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtFitter.h) + ⎿  #ifndef AtFITTER_H + #define AtFITTER_H + … +40 lines (ctrl+o to expand) + +● Bash(git show RealAurio/fitterRefactor:AtData/AtFittedTrack.h) + ⎿  #ifndef ATFITTEDTRACK_H + #define ATFITTEDTRACK_H + … +162 lines (ctrl+o to expand) + +● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitterTask.h && echo + "---" && git show + RealAurio/fitterRefactor:AtReconstruction/AtFitterTask.cxx) + ⎿  /********************************************************************* + * Fitter Task AtFitterTask.hh * + * Author: Y. Ayyad ayyadlim@frib.msu.edu * + … +95 lines (ctrl+o to expand) + +● Bash(git show RealAurio/fitterRefactor:AtData/AtFitTrackMetadata.h && echo + "---" && git show RealAurio/fitterRefactor:AtData/AtFitMetadata.h) + ⎿  #ifndef ATFITTRACKMETADATA_H + #define ATFITTRACKMETADATA_H + … +39 lines (ctrl+o to expand) + +● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtFitter.cxx) + ⎿  #include "AtFitter.h" + + #include "AtPatternEvent.h" + … +25 lines (ctrl+o to expand) + +● Searched for 2 patterns, read 1 file (ctrl+o to expand) + +● I now have the full picture. Here's the re-ask on output storage: + +● User answered Claude's questions: + ⎿  · What UKF output should be stored in AtFittedTrack beyond the + standard vertex state (Kinematics + Vertex + TrackProperties already + in the refactored class)? → Vertex + smoothed positions (Recommended) + +● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtGenfit.h | + head -60) + ⎿  #ifndef ATGENFIT_H + #define ATGENFIT_H + … +36 lines (ctrl+o to expand) + +● Bash(git show RealAurio/fitterRefactor:AtData/AtTrack.h | grep -n + "smoothed\|Smooth\|UKF\|fSmooth\|fFitted" | head -20) + ⎿  (No output) + +● Bash(git show + RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtFitterOld.h) + ⎿  #ifndef AtFITTEROLD_H + #define AtFITTEROLD_H + … +37 lines (ctrl+o to expand) + +● Bash(git show RealAurio/fitterRefactor:AtData/AtTrackingEventOld.h | head -30) + ⎿  #ifndef AtTRACKINGEVENTOLD_H + #define AtTRACKINGEVENTOLD_H + … +27 lines (ctrl+o to expand) + +● Searched for 1 pattern (ctrl+o to expand) + +● Now I have everything I need. Let me write the plan: + +● Updated plan + ⎿  /plan to preview + ⎿  User rejected Claude's plan: + ╭─────────────────────────────────────────────────────────────────────────╮ + │ Plan: Integrate UKF Track Fitter into Task System │ + │ │ + │ Context │ + │ │ + │ The UKF track fitter (kf::TrackFitterUKF) is fully implemented as a │ + │ standalone physics class and tested via macros in macro/tests/UKF/. It │ + │ needs to be integrated into the FairTask-based reconstruction pipeline │ + │ so it can be run event-by-event on data. │ + │ │ + │ In parallel, PR #254 (RealAurio/fitterRefactor) proposes a clean │ + │ refactored fitter architecture with a new EventFit::AtFitter base class │ + │ and a redesigned AtFittedTrack. The UKF task should be implemented │ + │ against this new interface rather than the deprecated │ + │ AtFITTER::AtFitter. │ + │ │ + │ Outcome: A complete EventFit::AtFitterUKF class that reads an │ + │ AtPatternEvent (via AtFitterTask), fits each track with the UKF using a │ + │ Brho-based momentum seed, and writes an AtTrackingEvent with │ + │ AtFittedTrack objects containing vertex kinematics and smoothed │ + │ trajectory positions. │ + │ │ + │ --- │ + │ Step 1 — Merge fitterRefactor (PR #254) into OpenKF-Claude │ + │ │ + │ Merge RealAurio/fitterRefactor into the OpenKF-Claude branch. This │ + │ brings in all 36 changed files. Key files introduced/changed: │ + │ │ + │ File: AtData/AtFittedTrack.h/cxx │ + │ Change: New nested-struct design (Kinematics, ParticleInfo, │ + │ TrackProperties) │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtFitMetadata.h/cxx │ + │ Change: New: event-level fit metadata │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtFitTrackMetadata.h/cxx │ + │ Change: New: per-track fit statistics │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtFittedTrackOld.h/cxx │ + │ Change: Backward-compat alias for old AtFittedTrack │ + │ ──────────────────────────────────────── │ + │ File: AtData/CMakeLists.txt │ + │ Change: Adds new data classes │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtDataLinkDef.h │ + │ Change: Adds ROOT dictionary entries │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtFitter.h/cxx │ + │ Change: New EventFit::AtFitter with FitEvent()/GetFittedTrack() │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtFitterOld.h/cxx │ + │ Change: Old AtFITTER::AtFitter kept for backward compat │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtGenfit.h │ + │ Change: Now inherits AtFitterOld (deprecated, not changed) │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtMCFitter.cxx │ + │ Change: Minor updates │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitterTask.h/cxx │ + │ Change: Updated to use EventFit::AtFitter; adds optional │ + │ RawEvent/Event/FitMetadata branches │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtReconstructionLinkDef.h │ + │ Change: Updated LinkDef │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtPatternModification/* │ + │ Change: New BraggCurveFinder task (independent) │ + │ ──────────────────────────────────────── │ + │ File: AtEventDisplay/AtTabs/* │ + │ Change: Tab updates (independent) │ + │ │ + │ Resolve any merge conflicts; most conflict risk is in CMakeLists and │ + │ LinkDef files. │ + │ │ + │ --- │ + │ Step 2 — Add Smoothed Positions to AtFittedTrack::TrackProperties │ + │ │ + │ File: AtData/AtFittedTrack.h │ + │ │ + │ Add to TrackProperties struct: │ + │ std::vector fSmoothedPositions; // UKF-fitted │ + │ positions at each cluster │ + │ │ + │ Add corresponding setter/getter to AtFittedTrack: │ + │ void SetSmoothedPositions(std::vector pos) { │ + │ fTrackProperties.fSmoothedPositions = std::move(pos); } │ + │ const std::vector& GetSmoothedPositions() const { │ + │ return fTrackProperties.fSmoothedPositions; } │ + │ │ + │ Bump ClassDef(AtFittedTrack, 3). │ + │ │ + │ --- │ + │ Step 3 — Implement AtFitterUKF │ + │ │ + │ New files: │ + │ - AtReconstruction/AtFitter/AtFitterUKF.h │ + │ - AtReconstruction/AtFitter/AtFitterUKF.cxx │ + │ │ + │ Header │ + │ │ + │ namespace EventFit { │ + │ │ + │ class AtFitterUKF : public AtFitter { │ + │ public: │ + │ AtFitterUKF(double charge, double mass_MeV, │ + │ std::unique_ptr elossModel); │ + │ │ + │ void SetBField(ROOT::Math::XYZVector bField); // Tesla │ + │ void SetEField(ROOT::Math::XYZVector eField); // V/m (default 0) │ + │ void SetUKFParameters(double alpha, double beta, double kappa); │ + │ void SetMeasurementSigma(double sigma_mm); // default 1.0 mm │ + │ void SetMomentumSigmaFrac(double frac); // fractional, default │ + │ 0.1 │ + │ void SetMinClusters(int n); // default 3 │ + │ │ + │ void Init() override {} │ + │ │ + │ protected: │ + │ AtFittedTrack *GetFittedTrack(AtTrack *track, │ + │ AtFitMetadata *fitMetadata = nullptr, │ + │ AtRawEvent *rawEvent = nullptr, │ + │ AtEvent *event = nullptr) override; │ + │ private: │ + │ // Helpers │ + │ std::unique_ptr CreateUKF() const; │ + │ ROOT::Math::XYZPoint GetInitialPosition(AtTrack *track) const; │ + │ ROOT::Math::XYZVector GetInitialMomentum(AtTrack *track) const; │ + │ TMatrixD GetInitialCovariance(double p_mag_MeV) const; │ + │ TMatrixD GetMeasCovariance() const; │ + │ double GetBrho(AtTrack *track) const; │ + │ │ + │ // Particle / detector config │ + │ double fCharge; // charge in │ + │ units of e │ + │ double fMass_MeV; // rest mass │ + │ in MeV/c² │ + │ ROOT::Math::XYZVector fBField{0, 0, 2.85}; // Tesla │ + │ ROOT::Math::XYZVector fEField{0, 0, 0}; // V/m │ + │ std::unique_ptr fELossModel; │ + │ │ + │ // UKF tuning │ + │ double fAlpha{1e-3}, fBeta{2.0}, fKappa{0.0}; │ + │ double fMeasSigma_mm{1.0}; │ + │ double fMomSigmaFrac{0.1}; │ + │ int fMinClusters{3}; │ + │ │ + │ AtTools::AtKinematics fKinematics; // for GetMomFromBrho │ + │ }; │ + │ │ + │ } // namespace EventFit │ + │ │ + │ Implementation logic in GetFittedTrack() │ + │ │ + │ 1. Get clusters: auto *clusters = track->GetHitClusterArray(); │ + │ track->SortClusterHitArrayZ(); │ + │ Skip track if clusters->size() < fMinClusters. │ + │ 2. Seed momentum from circle fit (same as GenFit): │ + │ radius = track->GetGeoRadius() / 1000.0 // mm → m │ + │ theta = track->GetGeoTheta() │ + │ phi = track->GetGeoPhi() │ + │ brho = fBField.Z() * radius / sin(theta) // T·m │ + │ [p_GeV, KE_GeV] = fKinematics.GetMomFromBrho(mass_amu, charge, brho) │ + │ p_MeV = p_GeV * 1000.0 │ + │ 3. Initial position: clusters->front().GetPosition() (in mm). │ + │ 4. Initial momentum direction: unit vector from (theta, phi); magnitude │ + │ = p_MeV. │ + │ 5. Create and configure UKF: │ + │ auto ukf = CreateUKF(); // new kf::TrackFitterUKF with │ + │ AtPropagator+AtRK4AdaptiveStepper │ + │ ukf->setParameters(fAlpha, fBeta, fKappa); │ + │ ukf->SetInitialState(initialPos, initialMom, │ + │ GetInitialCovariance(p_MeV)); │ + │ ukf->SetMeasCov(GetMeasCovariance()); │ + │ ukf->fEnableEnStraggling = true; │ + │ 6. Filter loop: │ + │ for (size_t i = 1; i < clusters->size(); ++i) { │ + │ XYZPoint measPt = clusters->at(i).GetPosition(); // mm │ + │ ukf->predictUKF(measPt); │ + │ ukf->correctUKF(measPt); │ + │ } │ + │ ukf->smoothUKF(); │ + │ 6. Wrap in try/catch for Cholesky failures → mark fitConverged = false. │ + │ 7. Build AtFittedTrack: │ + │ - From smoothedStates[0] = [x, y, z, p, theta_s, phi_s]: │ + │ - KE = sqrt(p*p + m*m) - m │ + │ - SetKinematics(KE, theta_s, phi_s) │ + │ - SetVertex(XYZVector{x, y, z}) │ + │ - SetParticleInfo(pdg_string, charge, mass_amu) │ + │ - Extract smoothed positions from all states → │ + │ SetSmoothedPositions(...) │ + │ - Compute trackLength as sum of distances between consecutive │ + │ smoothed positions │ + │ - SetTrackPropertiesStruct(...) with initialPosition, trackLength, │ + │ trackPoints = clusters->size() │ + │ - Create AtFitTrackMetadata: set trackID, fitConverged; chi² = sum of │ + │ squared smoothed residuals / ndf; SetTrackMetadata() │ + │ 8. Optionally store metadata in AtFitMetadata (if non-null) via │ + │ SetTrackMetadatasVector(trackID, ...). │ + │ │ + │ --- │ + │ Step 4 — Update Build System │ + │ │ + │ File: AtReconstruction/CMakeLists.txt │ + │ - Add AtFitter/AtFitterUKF.cxx to SRCS │ + │ │ + │ File: AtReconstruction/AtReconstructionLinkDef.h │ + │ - Add #pragma link C++ class EventFit::AtFitterUKF+; │ + │ │ + │ --- │ + │ Step 5 — Write Usage Macro │ + │ │ + │ New file: macro/tests/UKF/UKFTask.C │ + │ │ + │ Demonstrates: │ + │ auto eloss = std::make_unique(3.553e-5); │ + │ eloss->SetProjectile(1, 1, 1); │ + │ eloss->SetMaterial({{1, 1, 1}}); │ + │ │ + │ auto fitter = std::make_unique(1.0, 938.272, │ + │ std::move(eloss)); │ + │ fitter->SetBField({0, 0, 2.85}); │ + │ fitter->SetUKFParameters(1e-3, 2, 0); │ + │ fitter->SetMeasurementSigma(1.0); │ + │ │ + │ auto task = new AtFitterTask(std::move(fitter)); │ + │ task->SetPersistence(kTRUE); │ + │ run->AddTask(task); │ + │ │ + │ --- │ + │ MCFitter Assessment │ + │ │ + │ AtMCFitter uses stochastic multi-iteration Monte Carlo optimization and │ + │ produces AtMCResult arrays — not per-track AtFittedTrack objects. Its │ + │ Exec(AtPatternEvent&) pattern is fundamentally different from │ + │ GetFittedTrack(AtTrack*). It does not combine well with the new │ + │ EventFit::AtFitter base and should remain separate (AtMCFitterTask). │ + │ │ + │ --- │ + │ Verification │ + │ │ + │ 1. Build: cmake --build build -j$(nproc) — must compile without errors. │ + │ 2. Existing macro: Run macro/tests/UKF/UKFSingleTrack.C to confirm │ + │ TrackFitterUKF still works standalone. │ + │ 3. New macro: Run macro/tests/UKF/UKFTask.C end-to-end on a simulation │ + │ file; inspect the AtTrackingEvent branch in the output ROOT file — │ + │ verify AtFittedTrack contains non-trivial kinematics and smoothed │ + │ positions. │ + │ 4. Check merge: Verify AtGenfit still compiles (it uses AtFitterOld │ + │ which is unchanged). │ + ╰─────────────────────────────────────────────────────────────────────────╯ + +❯ To this plan, I think we should emphasize the need for proper tests using + the gtest stuff we use other places in the code-base. + +● Explore(Explore gtest structure in codebase) + ⎿  Done (30 tool uses · 52.4k tokens · 1m 4s) + (ctrl+o to expand) + +● Good. Let me also read the existing UKF test file to understand what's already + covered: + +● Read 1 file (ctrl+o to expand) + +● Now I have a complete picture of the test infrastructure. Let me update the + plan with a proper testing section: + +● Updated plan + ⎿  /plan to preview + ⎿  User rejected Claude's plan: + ╭─────────────────────────────────────────────────────────────────────────╮ + │ Plan: Integrate UKF Track Fitter into Task System │ + │ │ + │ Context │ + │ │ + │ The UKF track fitter (kf::TrackFitterUKF) is fully implemented as a │ + │ standalone physics class and tested via macros in macro/tests/UKF/. It │ + │ needs to be integrated into the FairTask-based reconstruction pipeline │ + │ so it can be run event-by-event on data. │ + │ │ + │ In parallel, PR #254 (RealAurio/fitterRefactor) proposes a clean │ + │ refactored fitter architecture with a new EventFit::AtFitter base class │ + │ and a redesigned AtFittedTrack. The UKF task should be implemented │ + │ against this new interface rather than the deprecated │ + │ AtFITTER::AtFitter. │ + │ │ + │ Outcome: A complete EventFit::AtFitterUKF class that reads an │ + │ AtPatternEvent (via AtFitterTask), fits each track with the UKF using a │ + │ Brho-based momentum seed, and writes an AtTrackingEvent with │ + │ AtFittedTrack objects containing vertex kinematics and smoothed │ + │ trajectory positions. │ + │ │ + │ --- │ + │ Step 1 — Merge fitterRefactor (PR #254) into OpenKF-Claude │ + │ │ + │ Merge RealAurio/fitterRefactor into the OpenKF-Claude branch. This │ + │ brings in all 36 changed files. Key files introduced/changed: │ + │ │ + │ File: AtData/AtFittedTrack.h/cxx │ + │ Change: New nested-struct design (Kinematics, ParticleInfo, │ + │ TrackProperties) │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtFitMetadata.h/cxx │ + │ Change: New: event-level fit metadata │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtFitTrackMetadata.h/cxx │ + │ Change: New: per-track fit statistics │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtFittedTrackOld.h/cxx │ + │ Change: Backward-compat alias for old AtFittedTrack │ + │ ──────────────────────────────────────── │ + │ File: AtData/CMakeLists.txt │ + │ Change: Adds new data classes │ + │ ──────────────────────────────────────── │ + │ File: AtData/AtDataLinkDef.h │ + │ Change: Adds ROOT dictionary entries │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtFitter.h/cxx │ + │ Change: New EventFit::AtFitter with FitEvent()/GetFittedTrack() │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtFitterOld.h/cxx │ + │ Change: Old AtFITTER::AtFitter kept for backward compat │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtGenfit.h │ + │ Change: Now inherits AtFitterOld (deprecated, not changed) │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitter/AtMCFitter.cxx │ + │ Change: Minor updates │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtFitterTask.h/cxx │ + │ Change: Updated to use EventFit::AtFitter; adds optional │ + │ RawEvent/Event/FitMetadata branches │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtReconstructionLinkDef.h │ + │ Change: Updated LinkDef │ + │ ──────────────────────────────────────── │ + │ File: AtReconstruction/AtPatternModification/* │ + │ Change: New BraggCurveFinder task (independent) │ + │ ──────────────────────────────────────── │ + │ File: AtEventDisplay/AtTabs/* │ + │ Change: Tab updates (independent) │ + │ │ + │ Resolve any merge conflicts; most conflict risk is in CMakeLists and │ + │ LinkDef files. │ + │ │ + │ --- │ + │ Step 2 — Add Smoothed Positions to AtFittedTrack::TrackProperties │ + │ │ + │ File: AtData/AtFittedTrack.h │ + │ │ + │ Add to TrackProperties struct: │ + │ std::vector fSmoothedPositions; // UKF-fitted │ + │ positions at each cluster │ + │ │ + │ Add corresponding setter/getter to AtFittedTrack: │ + │ void SetSmoothedPositions(std::vector pos) { │ + │ fTrackProperties.fSmoothedPositions = std::move(pos); } │ + │ const std::vector& GetSmoothedPositions() const { │ + │ return fTrackProperties.fSmoothedPositions; } │ + │ │ + │ Bump ClassDef(AtFittedTrack, 3). │ + │ │ + │ --- │ + │ Step 3 — Implement AtFitterUKF │ + │ │ + │ New files: │ + │ - AtReconstruction/AtFitter/AtFitterUKF.h │ + │ - AtReconstruction/AtFitter/AtFitterUKF.cxx │ + │ │ + │ Header │ + │ │ + │ namespace EventFit { │ + │ │ + │ class AtFitterUKF : public AtFitter { │ + │ public: │ + │ AtFitterUKF(double charge, double mass_MeV, │ + │ std::unique_ptr elossModel); │ + │ │ + │ void SetBField(ROOT::Math::XYZVector bField); // Tesla │ + │ void SetEField(ROOT::Math::XYZVector eField); // V/m (default 0) │ + │ void SetUKFParameters(double alpha, double beta, double kappa); │ + │ void SetMeasurementSigma(double sigma_mm); // default 1.0 mm │ + │ void SetMomentumSigmaFrac(double frac); // fractional, default │ + │ 0.1 │ + │ void SetMinClusters(int n); // default 3 │ + │ │ + │ void Init() override {} │ + │ │ + │ protected: │ + │ AtFittedTrack *GetFittedTrack(AtTrack *track, │ + │ AtFitMetadata *fitMetadata = nullptr, │ + │ AtRawEvent *rawEvent = nullptr, │ + │ AtEvent *event = nullptr) override; │ + │ private: │ + │ // Helpers │ + │ std::unique_ptr CreateUKF() const; │ + │ ROOT::Math::XYZPoint GetInitialPosition(AtTrack *track) const; │ + │ ROOT::Math::XYZVector GetInitialMomentum(AtTrack *track) const; │ + │ TMatrixD GetInitialCovariance(double p_mag_MeV) const; │ + │ TMatrixD GetMeasCovariance() const; │ + │ double GetBrho(AtTrack *track) const; │ + │ │ + │ // Particle / detector config │ + │ double fCharge; // charge in │ + │ units of e │ + │ double fMass_MeV; // rest mass │ + │ in MeV/c² │ + │ ROOT::Math::XYZVector fBField{0, 0, 2.85}; // Tesla │ + │ ROOT::Math::XYZVector fEField{0, 0, 0}; // V/m │ + │ std::unique_ptr fELossModel; │ + │ │ + │ // UKF tuning │ + │ double fAlpha{1e-3}, fBeta{2.0}, fKappa{0.0}; │ + │ double fMeasSigma_mm{1.0}; │ + │ double fMomSigmaFrac{0.1}; │ + │ int fMinClusters{3}; │ + │ │ + │ AtTools::AtKinematics fKinematics; // for GetMomFromBrho │ + │ }; │ + │ │ + │ } // namespace EventFit │ + │ │ + │ Implementation logic in GetFittedTrack() │ + │ │ + │ 1. Get clusters: auto *clusters = track->GetHitClusterArray(); │ + │ track->SortClusterHitArrayZ(); │ + │ Skip track if clusters->size() < fMinClusters. │ + │ 2. Seed momentum from circle fit (same as GenFit): │ + │ radius = track->GetGeoRadius() / 1000.0 // mm → m │ + │ theta = track->GetGeoTheta() │ + │ phi = track->GetGeoPhi() │ + │ brho = fBField.Z() * radius / sin(theta) // T·m │ + │ [p_GeV, KE_GeV] = fKinematics.GetMomFromBrho(mass_amu, charge, brho) │ + │ p_MeV = p_GeV * 1000.0 │ + │ 3. Initial position: clusters->front().GetPosition() (in mm). │ + │ 4. Initial momentum direction: unit vector from (theta, phi); magnitude │ + │ = p_MeV. │ + │ 5. Create and configure UKF: │ + │ auto ukf = CreateUKF(); // new kf::TrackFitterUKF with │ + │ AtPropagator+AtRK4AdaptiveStepper │ + │ ukf->setParameters(fAlpha, fBeta, fKappa); │ + │ ukf->SetInitialState(initialPos, initialMom, │ + │ GetInitialCovariance(p_MeV)); │ + │ ukf->SetMeasCov(GetMeasCovariance()); │ + │ ukf->fEnableEnStraggling = true; │ + │ 6. Filter loop: │ + │ for (size_t i = 1; i < clusters->size(); ++i) { │ + │ XYZPoint measPt = clusters->at(i).GetPosition(); // mm │ + │ ukf->predictUKF(measPt); │ + │ ukf->correctUKF(measPt); │ + │ } │ + │ ukf->smoothUKF(); │ + │ 6. Wrap in try/catch for Cholesky failures → mark fitConverged = false. │ + │ 7. Build AtFittedTrack: │ + │ - From smoothedStates[0] = [x, y, z, p, theta_s, phi_s]: │ + │ - KE = sqrt(p*p + m*m) - m │ + │ - SetKinematics(KE, theta_s, phi_s) │ + │ - SetVertex(XYZVector{x, y, z}) │ + │ - SetParticleInfo(pdg_string, charge, mass_amu) │ + │ - Extract smoothed positions from all states → │ + │ SetSmoothedPositions(...) │ + │ - Compute trackLength as sum of distances between consecutive │ + │ smoothed positions │ + │ - SetTrackPropertiesStruct(...) with initialPosition, trackLength, │ + │ trackPoints = clusters->size() │ + │ - Create AtFitTrackMetadata: set trackID, fitConverged; chi² = sum of │ + │ squared smoothed residuals / ndf; SetTrackMetadata() │ + │ 8. Optionally store metadata in AtFitMetadata (if non-null) via │ + │ SetTrackMetadatasVector(trackID, ...). │ + │ │ + │ --- │ + │ Step 4 — Update Build System │ + │ │ + │ File: AtReconstruction/CMakeLists.txt │ + │ - Add AtFitter/AtFitterUKF.cxx to SRCS │ + │ │ + │ File: AtReconstruction/AtReconstructionLinkDef.h │ + │ - Add #pragma link C++ class EventFit::AtFitterUKF+; │ + │ │ + │ --- │ + │ Step 5 — GTest Unit Tests │ + │ │ + │ New file: AtReconstruction/AtFitter/AtFitterUKFTest.cxx │ + │ │ + │ Register in AtReconstruction/CMakeLists.txt inside the existing │ + │ if(TARGET Eigen3::Eigen) block (same block as TrackFitterUKFTest.cxx): │ + │ │ + │ set(OPENK_TEST_SRCS │ + │ ... │ + │ AtFitter/TrackFitterUKFTest.cxx # already there │ + │ AtFitter/AtFitterUKFTest.cxx # new │ + │ ) │ + │ attpcroot_generate_tests(OpenKFTests SRCS ${OPENK_TEST_SRCS} DEPS │ + │ AtReconstruction) │ + │ │ + │ Use the same test-class subclassing pattern as TrackFitterUKFTest.cxx: │ + │ define AtFitterUKFTestable that exposes internals. │ + │ │ + │ Test fixture │ + │ │ + │ class AtFitterUKFFixture : public testing::Test { │ + │ protected: │ + │ // Known proton track from simulation (same data used in │ + │ UKFSingleTrack.C macro) │ + │ // Start position in mm, momentum in MeV/c │ + │ static constexpr double mass_p = 938.272; │ + │ static constexpr double charge_p = 1.602176634e-19; // C │ + │ static constexpr double B_z = 2.85; // T │ + │ static constexpr double density = 3.553e-5; // g/cm³ H2 │ + │ │ + │ std::unique_ptr fitter; │ + │ std::vector fMeasurements; // clusters from │ + │ simulation │ + │ │ + │ void SetUp() override { │ + │ // Build ELoss model │ + │ auto eloss = std::make_unique(0); │ + │ eloss->LoadSrimTable(getEnergyPath()); │ + │ │ + │ fitter = std::make_unique(charge_p, │ + │ mass_p, std::move(eloss)); │ + │ fitter->SetBField({0, 0, B_z}); │ + │ fitter->SetUKFParameters(1e-3, 2, 0); │ + │ fitter->SetMeasurementSigma(1.0); │ + │ │ + │ // Load known cluster positions from hits.txt (in mm) │ + │ fMeasurements = LoadHitsFromFile(getHitsPath()); │ + │ } │ + │ │ + │ // Build a synthetic AtTrack from the known measurements │ + │ AtTrack BuildTrack(); │ + │ }; │ + │ │ + │ Tests to implement │ + │ │ + │ ┌─────────────────────────────────────┬──────────────────────────────── │ + │ ───┐ │ + │ │ Test name │ What is checked │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ Construction │ AtFitterUKF constructs without │ + │ │ │ + │ │ │ throwing │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ SkipsTrackWithFewClusters │ Returns nullptr when cluster │ + │ │ │ + │ │ │ count < fMinClusters │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ │ Fits synthetic proton track; │ + │ │ │ + │ │ FittedTrackHasReasonableKinematics │ asserts KE > 0, 0 < theta < π, │ + │ │ │ + │ │ │ smoothed positions non-empty │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ FittedTrackKinematicsWithinToleranc │ Asserts fitted KE is within 20% │ + │ │ │ + │ │ e │ of known true value (proton at │ + │ │ │ + │ │ │ ~47 MeV/c from simulation) │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ SmoothedPositionsPopulated │ GetSmoothedPositions().size() │ + │ == │ │ + │ │ │ nClusters - 1 │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ MetadataConvergedSet │ │ + │ GetTrackMetadata()->GetFitConverg │ │ + │ │ │ ed() == true for good track │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ │ SetKinematics(E,θ,φ) → │ + │ │ │ + │ │ AtFittedTrackKinematicsRoundTrip │ GetKinematics() returns the │ + │ same │ │ + │ │ │ values │ + │ │ │ + │ ├─────────────────────────────────────┼──────────────────────────────── │ + │ ───┤ │ + │ │ AtFittedTrackSmoothedPositionsRound │ SetSmoothedPositions(v) → │ + │ │ │ + │ │ Trip │ GetSmoothedPositions() returns │ + │ │ │ + │ │ │ same vector │ + │ │ │ + │ └─────────────────────────────────────┴──────────────────────────────── │ + │ ───┘ │ + │ │ + │ The hits.txt file already exists in macro/tests/UKF/hits.txt and │ + │ contains the simulation cluster positions. Use $VMCWORKDIR env var or a │ + │ relative fallback path (same pattern as getEnergyPath() in │ + │ TrackFitterUKFTest.cxx). │ + │ │ + │ --- │ + │ Step 6 — Write Usage Macro │ + │ │ + │ New file: macro/tests/UKF/UKFTask.C │ + │ │ + │ Demonstrates: │ + │ auto eloss = std::make_unique(3.553e-5); │ + │ eloss->SetProjectile(1, 1, 1); │ + │ eloss->SetMaterial({{1, 1, 1}}); │ + │ │ + │ auto fitter = std::make_unique(charge_p, mass_p, │ + │ std::move(eloss)); │ + │ fitter->SetBField({0, 0, 2.85}); │ + │ fitter->SetUKFParameters(1e-3, 2, 0); │ + │ fitter->SetMeasurementSigma(1.0); │ + │ │ + │ auto task = new AtFitterTask(std::move(fitter)); │ + │ task->SetPersistence(kTRUE); │ + │ run->AddTask(task); │ + │ │ + │ --- │ + │ MCFitter Assessment │ + │ │ + │ AtMCFitter uses stochastic multi-iteration Monte Carlo optimization and │ + │ produces AtMCResult arrays — not per-track AtFittedTrack objects. Its │ + │ Exec(AtPatternEvent&) pattern is fundamentally different from │ + │ GetFittedTrack(AtTrack*). It does not combine well with the new │ + │ EventFit::AtFitter base and should remain separate (AtMCFitterTask). │ + │ │ + │ --- │ + │ Verification │ + │ │ + │ 1. Build: cmake --build build -j$(nproc) — must compile without errors. │ + │ 2. Existing macro: Run macro/tests/UKF/UKFSingleTrack.C to confirm │ + │ TrackFitterUKF still works standalone. │ + │ 3. New macro: Run macro/tests/UKF/UKFTask.C end-to-end on a simulation │ + │ file; inspect the AtTrackingEvent branch in the output ROOT file — │ + │ verify AtFittedTrack contains non-trivial kinematics and smoothed │ + │ positions. │ + │ 4. Check merge: Verify AtGenfit still compiles (it uses AtFitterOld │ + │ which is unchanged). │ + ╰─────────────────────────────────────────────────────────────────────────╯ + +❯ I just added your CLAUDE.md file. I didn't realize it was missing. Please + update your plan with the knowledge that you will always have that + context + +● Searched for 1 pattern, read 3 files (ctrl+o to expand) + +● Key impacts from CLAUDE.md on the plan: + 1. Tests must not access external files — the hits.txt loading approach is + invalid; need to hardcode cluster data + 2. LinkDef suffix rules — +; for persisted classes, -!; for tasks/algorithms + 3. C++ types in non-persisted classes — use double/int/bool, not + Double_t/Int_t in AtFitterUKF + 4. Build command is cmake --build build -j10 + 5. clang-format-17 must be run on all new files + +● Updated plan + ⎿  /plan to preview + +● Updated plan + ⎿  /plan to preview + +● Now fix the test section to remove file access and fix the build command in + verification: + +● Updated plan + ⎿  /plan to preview + +● Now fix the verification section to use the correct build command: + +● Updated plan + ⎿  /plan to preview + +● Also save relevant CLAUDE.md facts to memory before exiting plan mode: + +● Recalled 1 memory, wrote 1 memory (ctrl+o to expand) + ⎿  Tool use rejected with user message: /export + +● Skill(export) + ⎿  Initializing… + ⎿  Error: Skill export is not a prompt-based skill + +● I'm not familiar with a /export command in this context. Could you clarify + what you'd like to export, or what you meant? + +✻ Crunched for 3m 10s \ No newline at end of file diff --git a/.claude/CLAUDE.md b/.claude/CLAUDE.md index 6e3f1de01..138117faa 100644 --- a/.claude/CLAUDE.md +++ b/.claude/CLAUDE.md @@ -1,146 +1,60 @@ # CLAUDE.md -This file provides guidance to Claude Code (claude.ai/code) when working with code in this repository. - -## About - -ATTPCROOT is a ROOT/FairRoot-based C++ framework for simulation and analysis of Active Target Time Projection Chamber (AT-TPC) detector data. It integrates with FairSoft (provides ROOT, Geant4, VMC) and FairRoot (provides the task/run framework). - -## Environment Setup - -Before building or running anything, the environment must be loaded. The VSCode terminal auto-sources this on startup: - -```bash -source build/config.sh -``` - -This sets `LD_LIBRARY_PATH`, `ROOTSYS`, `VMCWORKDIR`, `ROOT_INCLUDE_PATH`, and Geant4 data paths. The key install paths on this machine are: -- FairRoot: `~/fair_install/FairRootInstall` -- FairSoft: `~/fair_install/FairSoftInstall` - -For CMake configuration (not the build itself), the following env vars are needed: -```bash -export FAIRROOTPATH=~/fair_install/FairRootInstall -export SIMPATH=~/fair_install/FairSoftInstall -``` - -## Build Commands - -```bash -# Configure (from repo root, out-of-source into build/) -cmake -S . -B build -DCMAKE_PREFIX_PATH=~/fair_install/hdf5 - -# Build (10 parallel jobs) -cmake --build build -j10 - -# Build a specific target -cmake --build build --target AtSimulationData -j10 - -# Run all unit tests -cd build && ctest -V - -# Run a specific test binary directly -./build/tests/AtSimulationDataTests -./build/tests/AtGeneratorsTests -./build/tests/AtToolsTests -``` - -Tests are built by default (`BUILD_TESTS=ON`). Test binaries are placed in `build/tests/`. - -## Code Formatting - -The project uses `clang-format-17` with the config in `.clang-format` (based on LLVM style, 3-space indent, 120-char column limit). Format on save is configured in VSCode. To format manually: +ATTPCROOT is a ROOT/FairRoot-based C++ framework for simulation and analysis of Active Target Time Projection Chamber (AT-TPC) detector data. + +## Documentation + +Full developer documentation lives in `docs/`. See [docs/index.md](../docs/index.md) for the full map. Quick topic links: + +| Topic | File | +|-------|------| +| First-time install | [tooling/installation.md](../docs/tooling/installation.md) | +| Daily use (build/test) | [tooling/daily-use.md](../docs/tooling/daily-use.md) | +| Testing patterns | [tooling/testing.md](../docs/tooling/testing.md) | +| Contributor guide | [contributing/guide.md](../docs/contributing/guide.md) | +| Adding a new module | [contributing/new-module.md](../docs/contributing/new-module.md) | +| Code style | [contributing/code-style.md](../docs/contributing/code-style.md) | +| Module overview | [reference/modules.md](../docs/reference/modules.md) | +| Data model | [reference/data-model.md](../docs/reference/data-model.md) | +| Branch I/O contracts | [reference/branch-io-contracts.md](../docs/reference/branch-io-contracts.md) | +| Simulation pipeline | [subsystems/simulation-pipeline.md](../docs/subsystems/simulation-pipeline.md) | +| Reconstruction pipeline | [subsystems/reconstruction-pipeline.md](../docs/subsystems/reconstruction-pipeline.md) | +| Event generators | [subsystems/generators.md](../docs/subsystems/generators.md) | +| Pulse shape analysis | [subsystems/psa.md](../docs/subsystems/psa.md) | +| Energy loss | [subsystems/energy-loss.md](../docs/subsystems/energy-loss.md) | + +## Quick Reference: Build & Test ```bash -clang-format-17 -i -# Or format all changed files: -scripts/formatAll.sh +source build/config.sh # load environment (do this first) +cmake --build build -j10 # build everything +cd build && ctest -V # run all unit tests ``` -Static analysis uses `clang-tidy` (config in `.clang-tidy`): modernize-* and cppcoreguidelines-* checks. - -## Architecture - -The framework follows FairRoot's task pipeline pattern: data flows through a chain of `FairTask` subclasses, persisted as ROOT TClonesArrays in a TTree. - -### Module Overview - -| Module | Purpose | -|--------|---------| -| `AtSimulationData` | MC truth data: `AtStack`, `AtMCTrack`, `AtMCPoint`, `AtVertexPropagator` | -| `AtGenerators` | FairGenerator subclasses for beam/reaction/decay event generation | -| `AtDetectors` | Geant4 sensitive detector implementations (AT-TPC, GADGET-II, SpecMAT, etc.) | -| `AtDigitization` | Converts MC points → simulated pad signals (`AtClusterize`, `AtPulse`) | -| `AtUnpack` | Unpacks raw experimental data (GET/GRAW, HDF5, ROOT formats) | -| `AtData` | Core data classes: `AtRawEvent`, `AtEvent`, `AtHit`, `AtPad`, `AtTrack` | -| `AtMap` | Pad mapping between electronics channels and detector geometry | -| `AtParameter` | FairRuntimeDb parameter containers | -| `AtReconstruction` | PSA, pattern recognition, fitting tasks | -| `AtTools` | Utilities: energy loss (CATIMA), kinematics, space charge, hit sampling | -| `AtAnalysis` | High-level analysis tasks and `AtRunAna` | -| `AtS800` | S800 spectrograph data handling | -| `AtEventDisplay` | ROOT Eve-based event display | - -### Data Flow - -**Simulation pipeline** (macros in `macro/Simulation/`): -1. `FairPrimaryGenerator` with an `AtReactionGenerator` subclass → particle stack -2. Geant4/VMC transport → `AtMCPoint` hits in sensitive volumes -3. `AtClusterizeTask` → electron clusters from ionization -4. `AtPulseTask` → simulated pad signals (`AtRawEvent`) - -**Reconstruction pipeline** (macros in `macro/Unpack_*/` and `macro/Analysis/`): -1. `AtUnpackTask` → `AtRawEvent` (raw pad traces) -2. `AtFilterTask` → filtered `AtRawEvent` -3. `AtPSAtask` (Pulse Shape Analysis) → `AtEvent` with `AtHit` objects -4. `AtSampleConsensusTask` / `AtPRAtask` → `AtPatternEvent` with `AtTrack` objects -5. `AtMCFitterTask` / `AtFitterTask` → fitted tracks with kinematics - -### Key Design Patterns - -**FairTask subclasses**: Each processing step is a `FairTask`. They retrieve input branches via `TClonesArray*` from `FairRootManager` in `Init()` and process them in `Exec()`. - -**AtReactionGenerator**: Abstract base for all reaction generators. Subclasses implement `GenerateReaction()`. The `ReadEvent()` method is `final` and handles the beam/reaction event alternation via `AtVertexPropagator`. Generators can be chained for sequential decays using `SetSequentialDecay(true)`. - -**AtVertexPropagator**: Singleton (`AtVertexPropagator::Instance()`) that communicates vertex, momentum, and kinematics between chained generators and downstream tasks. Alternates beam/reaction events via `EndEvent()`. Use `ResetForTesting()` in unit tests. - -**PSA (Pulse Shape Analysis)**: `AtPSA` is the abstract base. Concrete implementations (`AtPSAMax`, `AtPSAFull`, `AtPSADeconv`, etc.) extract hits from pad traces. `AtPSAComposite` allows chaining PSA methods. - -**Energy loss**: `AtELossModel` is the base; `AtELossCATIMA` wraps the CATIMA library (fetched automatically if not installed). Used via `AtELossManager`. - -### Adding a New Module Library - -Each module's `CMakeLists.txt` follows this pattern: -```cmake -set(LIBRARY_NAME MyModule) -set(SRCS file1.cxx file2.cxx) -set(DEPENDENCIES ATTPCROOT::AtData ROOT::Core ...) -set(TEST_SRCS MyModuleTest.cxx) -attpcroot_generate_tests(${LIBRARY_NAME}Tests SRCS ${TEST_SRCS} DEPS ${LIBRARY_NAME}) -generate_target_and_root_library(${LIBRARY_NAME} LINKDEF ${LIBRARY_NAME}LinkDef.h SRCS ${SRCS} DEPS_PUBLIC ${DEPENDENCIES}) -``` +## Code-Writing Rules -ROOT dictionary generation requires a `*LinkDef.h` file listing every class that needs ROOT reflection. Every class with `ClassDef` in its header must appear in the LinkDef. +These rules apply whenever editing or adding C++ code in this repo. -**LinkDef streamer suffix** — the suffix after the class name controls whether ROOT generates a full I/O streamer: +### LinkDef Streamer Suffixes -- `ClassName +;` — generates a full streamer. Use this **only** for classes that are written to a ROOT file (i.e. stored as a branch in a `TClonesArray` or `TTree`). Examples: `AtHit`, `AtEvent`, `AtRawEvent`, `AtMCPoint`. -- `ClassName -!;` — registers the class for reflection (usable in interpreted macros, usable as a pointer type) but **suppresses the streamer**. Use this for every class that is never written to disk: tasks, algorithms, models, samplers, etc. Examples: `AtELossModel`, `AtSpaceChargeModel`, `AtSample` subclasses. +Every class in a `*LinkDef.h` file must use the correct suffix: -Generating an unnecessary streamer bloats the dictionary and binary with unused I/O code, so the default should be `-!` unless disk persistence is actually required. +- `ClassName +;` — generates a full I/O streamer. **Only** for classes written to a ROOT file (stored in a `TClonesArray` or `TTree` branch). Examples: `AtHit`, `AtEvent`, `AtRawEvent`, `AtMCPoint`. +- `ClassName -!;` — reflection only, no streamer. Use for **everything else**: tasks, algorithms, models, samplers. Examples: `AtPSAMax`, `AtELossModel`, `AtFitterTask`. -**C++ vs ROOT typedefs** — use plain C++ types (`bool`, `int`, `double`, `std::string`) in all classes that are not written to a ROOT file. Only use ROOT's portability typedefs (`Bool_t`, `Int_t`, `Double_t`, etc.) in classes that are persisted to disk, where ROOT's cross-platform I/O guarantees matter. Mixing ROOT types into purely algorithmic or task classes is unnecessary overhead. +Default to `-!` unless disk persistence is actually required. -### Unit Tests +### C++ vs ROOT Typedefs -Tests use Google Test, placed alongside source files (e.g., `AtVertexPropagatorTest.cxx` in `AtSimulationData/`). Register them in the module's `CMakeLists.txt` via `attpcroot_generate_tests()`. Tests must not access external files or network resources. +- Non-persisted classes (tasks, algorithms, models): use `bool`, `int`, `double`, `std::string`. +- Persisted data classes only: use `Bool_t`, `Int_t`, `Double_t`, etc. -### Macros +### Test Isolation -ROOT macros (`.C` files) in `macro/` are the primary user interface for running simulations and analyses. They are not compiled into libraries—they run interpreted by ROOT. Integration tests live in `macro/tests/`. +Unit tests must not access external files or network resources. Hardcode test data inline. ## Contributing -- PRs must target the `develop` branch and apply cleanly (fast-forward, no merge commits). -- Commit messages: present imperative mood, ≤72 characters summary. -- All PRs are checked by `clang-format`, `clang-tidy`, and the unit test suite. +- PRs target the `develop` branch; fast-forward only (no merge commits). +- Commit messages: present imperative mood, ≤72 characters. +- All PRs must pass `clang-format`, `clang-tidy`, and unit tests. diff --git a/.claude/plans/keen-marinating-wind.md b/.claude/plans/keen-marinating-wind.md new file mode 100644 index 000000000..f5693c34a --- /dev/null +++ b/.claude/plans/keen-marinating-wind.md @@ -0,0 +1,295 @@ +# Plan: Integrate UKF Track Fitter into Task System + +## Context +The UKF track fitter (`kf::TrackFitterUKF`) is fully implemented as a standalone physics class and tested via macros in `macro/tests/UKF/`. It needs to be integrated into the FairTask-based reconstruction pipeline so it can be run event-by-event on data. + +In parallel, PR #254 (`RealAurio/fitterRefactor`) proposes a clean refactored fitter architecture with a new `EventFit::AtFitter` base class and a redesigned `AtFittedTrack`. The UKF task should be implemented against this new interface rather than the deprecated `AtFITTER::AtFitter`. + +**Outcome:** A complete `EventFit::AtFitterUKF` class that reads an `AtPatternEvent` (via `AtFitterTask`), fits each track with the UKF using a Brho-based momentum seed, and writes an `AtTrackingEvent` with `AtFittedTrack` objects containing vertex kinematics and smoothed trajectory positions. + +--- + +## Step 1 — Merge fitterRefactor (PR #254) into OpenKF-Claude + +Merge `RealAurio/fitterRefactor` into the `OpenKF-Claude` branch. This brings in all 36 changed files. Key files introduced/changed: + +| File | Change | +|------|--------| +| `AtData/AtFittedTrack.h/cxx` | New nested-struct design (`Kinematics`, `ParticleInfo`, `TrackProperties`) | +| `AtData/AtFitMetadata.h/cxx` | New: event-level fit metadata | +| `AtData/AtFitTrackMetadata.h/cxx` | New: per-track fit statistics | +| `AtData/AtFittedTrackOld.h/cxx` | Backward-compat alias for old AtFittedTrack | +| `AtData/CMakeLists.txt` | Adds new data classes | +| `AtData/AtDataLinkDef.h` | Adds ROOT dictionary entries | +| `AtReconstruction/AtFitter/AtFitter.h/cxx` | New `EventFit::AtFitter` with `FitEvent()`/`GetFittedTrack()` | +| `AtReconstruction/AtFitter/AtFitterOld.h/cxx` | Old `AtFITTER::AtFitter` kept for backward compat | +| `AtReconstruction/AtFitter/AtGenfit.h` | Now inherits `AtFitterOld` (deprecated, not changed) | +| `AtReconstruction/AtFitter/AtMCFitter.cxx` | Minor updates | +| `AtReconstruction/AtFitterTask.h/cxx` | Updated to use `EventFit::AtFitter`; adds optional RawEvent/Event/FitMetadata branches | +| `AtReconstruction/AtReconstructionLinkDef.h` | Updated LinkDef | +| `AtReconstruction/AtPatternModification/*` | New BraggCurveFinder task (independent) | +| `AtEventDisplay/AtTabs/*` | Tab updates (independent) | + +Resolve any merge conflicts; most conflict risk is in CMakeLists and LinkDef files. + +--- + +## Step 2 — Add Smoothed Positions to `AtFittedTrack::TrackProperties` + +**File:** `AtData/AtFittedTrack.h` + +Add to `TrackProperties` struct: +```cpp +std::vector fSmoothedPositions; // UKF-fitted positions at each cluster +``` + +Add corresponding setter/getter to `AtFittedTrack`: +```cpp +void SetSmoothedPositions(std::vector pos) { fTrackProperties.fSmoothedPositions = std::move(pos); } +const std::vector& GetSmoothedPositions() const { return fTrackProperties.fSmoothedPositions; } +``` + +Bump `ClassDef(AtFittedTrack, 3)`. + +--- + +## Step 3 — Implement `AtFitterUKF` + +**New files:** +- `AtReconstruction/AtFitter/AtFitterUKF.h` +- `AtReconstruction/AtFitter/AtFitterUKF.cxx` + +**Type convention (from CLAUDE.md):** `AtFitterUKF` is never written to disk, so use plain C++ types (`double`, `int`, `bool`) — not ROOT typedefs (`Double_t`, `Int_t`, `Bool_t`). + +### Header + +```cpp +namespace EventFit { + +class AtFitterUKF : public AtFitter { +public: + AtFitterUKF(double charge, double mass_MeV, + std::unique_ptr elossModel); + + void SetBField(ROOT::Math::XYZVector bField); // Tesla + void SetEField(ROOT::Math::XYZVector eField); // V/m (default 0) + void SetUKFParameters(double alpha, double beta, double kappa); + void SetMeasurementSigma(double sigma_mm); // default 1.0 mm + void SetMomentumSigmaFrac(double frac); // fractional, default 0.1 + void SetMinClusters(int n); // default 3 + + void Init() override {} + +protected: + AtFittedTrack *GetFittedTrack(AtTrack *track, + AtFitMetadata *fitMetadata = nullptr, + AtRawEvent *rawEvent = nullptr, + AtEvent *event = nullptr) override; +private: + // Helpers + std::unique_ptr CreateUKF() const; + ROOT::Math::XYZPoint GetInitialPosition(AtTrack *track) const; + ROOT::Math::XYZVector GetInitialMomentum(AtTrack *track) const; + TMatrixD GetInitialCovariance(double p_mag_MeV) const; + TMatrixD GetMeasCovariance() const; + double GetBrho(AtTrack *track) const; + + // Particle / detector config (plain C++ types — not persisted) + double fCharge; + double fMass_MeV; + ROOT::Math::XYZVector fBField{0, 0, 2.85}; + ROOT::Math::XYZVector fEField{0, 0, 0}; + std::unique_ptr fELossModel; + + // UKF tuning + double fAlpha{1e-3}, fBeta{2.0}, fKappa{0.0}; + double fMeasSigma_mm{1.0}; + double fMomSigmaFrac{0.1}; + int fMinClusters{3}; + + AtTools::AtKinematics fKinematics; // for GetMomFromBrho + + ClassDef(AtFitterUKF, 1); +}; + +} // namespace EventFit +``` + +### Implementation logic in `GetFittedTrack()` + +1. **Get clusters**: `auto *clusters = track->GetHitClusterArray(); track->SortClusterHitArrayZ();` + Skip track if `clusters->size() < fMinClusters`. + +2. **Seed momentum from circle fit** (same as GenFit): + ``` + radius = track->GetGeoRadius() / 1000.0 // mm → m + theta = track->GetGeoTheta() + phi = track->GetGeoPhi() + brho = fBField.Z() * radius / sin(theta) // T·m + [p_GeV, KE_GeV] = fKinematics.GetMomFromBrho(mass_amu, charge, brho) + p_MeV = p_GeV * 1000.0 + ``` + +3. **Initial position**: `clusters->front().GetPosition()` (in mm). + +4. **Initial momentum direction**: unit vector from (theta, phi); magnitude = `p_MeV`. + +5. **Create and configure UKF**: + ```cpp + auto ukf = CreateUKF(); // new kf::TrackFitterUKF with AtPropagator+AtRK4AdaptiveStepper + ukf->setParameters(fAlpha, fBeta, fKappa); + ukf->SetInitialState(initialPos, initialMom, GetInitialCovariance(p_MeV)); + ukf->SetMeasCov(GetMeasCovariance()); + ukf->fEnableEnStraggling = true; + ``` + +6. **Filter loop**: + ```cpp + for (size_t i = 1; i < clusters->size(); ++i) { + XYZPoint measPt = clusters->at(i).GetPosition(); // mm + ukf->predictUKF(measPt); + ukf->correctUKF(measPt); + } + ukf->smoothUKF(); + ``` + Wrap in try/catch for Cholesky failures → mark `fitConverged = false`. + +7. **Build `AtFittedTrack`**: + - From `smoothedStates[0] = [x, y, z, p, theta_s, phi_s]`: + - `KE = sqrt(p*p + m*m) - m` + - `SetKinematics(KE, theta_s, phi_s)` + - `SetVertex(XYZVector{x, y, z})` + - `SetParticleInfo(pdg_string, charge, mass_amu)` + - Extract smoothed positions from all states → `SetSmoothedPositions(...)` + - Compute `trackLength` as sum of distances between consecutive smoothed positions + - `SetTrackPropertiesStruct(...)` with `initialPosition`, `trackLength`, `trackPoints = clusters->size()` + - Create `AtFitTrackMetadata`: set `trackID`, `fitConverged`; chi² = sum of squared smoothed residuals / `ndf`; `SetTrackMetadata()` + +8. **Optionally store metadata** in `AtFitMetadata` (if non-null) via `SetTrackMetadatasVector(trackID, ...)`. + +--- + +## Step 4 — Update Build System and LinkDef + +**File:** `AtReconstruction/CMakeLists.txt` +- Add `AtFitter/AtFitterUKF.cxx` to `SRCS` + +**File:** `AtReconstruction/AtReconstructionLinkDef.h` +- Add `#pragma link C++ class EventFit::AtFitterUKF-!;` — not persisted to disk; suppress streamer + +**File:** `AtData/AtDataLinkDef.h` (already updated by Step 1 merge, but verify): +- `AtFittedTrack +;` — persisted, needs full streamer +- `AtFitMetadata +;`, `AtFitTrackMetadata +;` — persisted +- `AtFittedTrack::TrackProperties +;` — nested struct in a persisted class needs its own entry +- `AtFittedTrackOld -!;`, etc. — backward-compat wrappers, no streamer needed + +**Code formatting:** Run `clang-format-17 -i` on all new/modified `.h` and `.cxx` files before committing, or use `scripts/formatAll.sh`. + +--- + +## Step 5 — GTest Unit Tests + +**New file:** `AtReconstruction/AtFitter/AtFitterUKFTest.cxx` + +Register in `AtReconstruction/CMakeLists.txt` inside the existing `if(TARGET Eigen3::Eigen)` block (same block as `TrackFitterUKFTest.cxx`): + +```cmake +set(OPENK_TEST_SRCS + ... + AtFitter/TrackFitterUKFTest.cxx # already there + AtFitter/AtFitterUKFTest.cxx # new +) +attpcroot_generate_tests(OpenKFTests SRCS ${OPENK_TEST_SRCS} DEPS AtReconstruction) +``` + +Use the same test-class subclassing pattern as `TrackFitterUKFTest.cxx`: define `AtFitterUKFTestable` that exposes internals. + +### Test fixture + +From CLAUDE.md: **"Tests must not access external files or network resources."** Cluster positions are therefore hardcoded directly from the known simulation track (same values used in `UKFSingleTrack.C`). The energy loss model uses `AtELossTable` loaded via `$VMCWORKDIR` — the same `getEnergyPath()` pattern already present in `TrackFitterUKFTest.cxx` (file is part of the repo, not external). + +```cpp +class AtFitterUKFFixture : public testing::Test { +protected: + static constexpr double mass_p = 938.272; // MeV/c² + static constexpr double charge_p = 1.602176634e-19; // C + static constexpr double B_z = 2.85; // T + + std::unique_ptr fitter; + + // Cluster positions hardcoded from simulation (mm) — no file access + const std::vector kClusters = { + {-3.4e-4, -1.5e-4, 1.0018}, // converted from cm, first point + {-14.895, -48.787, 10.1217}, + // ... (5–8 representative points from the known proton track) + }; + + void SetUp() override { + auto eloss = std::make_unique(0); + eloss->LoadSrimTable(getEnergyPath()); // repo file via VMCWORKDIR + eloss->SetDensity(3.553e-5); + + fitter = std::make_unique(charge_p, mass_p, std::move(eloss)); + fitter->SetBField({0, 0, B_z}); + fitter->SetUKFParameters(1e-3, 2, 0); + fitter->SetMeasurementSigma(1.0); + } + + // Build AtTrack with hardcoded clusters + PRA-derived Brho seed + AtTrack BuildTrack(bool tooFewClusters = false); +}; +``` + +### Tests to implement + +| Test name | What is checked | +|-----------|----------------| +| `Construction` | `AtFitterUKF` constructs without throwing | +| `SkipsTrackWithFewClusters` | Returns `nullptr` when cluster count < `fMinClusters` | +| `FittedTrackHasReasonableKinematics` | Fits synthetic proton track; `KE > 0`, `0 < theta < π`, smoothed positions non-empty | +| `FittedTrackKinematicsWithinTolerance` | Fitted KE within 20% of known true value (~47 MeV total for simulation proton) | +| `SmoothedPositionsPopulated` | `GetSmoothedPositions().size() == nClusters - 1` | +| `MetadataConvergedSet` | `GetTrackMetadata()->GetFitConverged() == true` for good track | +| `AtFittedTrackKinematicsRoundTrip` | `SetKinematics(E,θ,φ)` → `GetKinematics()` returns same values (no physics needed) | +| `AtFittedTrackSmoothedPositionsRoundTrip` | `SetSmoothedPositions(v)` → `GetSmoothedPositions()` returns same vector (no physics needed) | + +The round-trip tests (`AtFittedTrack*`) require no energy loss model and test the data class in isolation. + +--- + +## Step 6 — Write Usage Macro + +**New file:** `macro/tests/UKF/UKFTask.C` + +Demonstrates: +```cpp +auto eloss = std::make_unique(3.553e-5); +eloss->SetProjectile(1, 1, 1); +eloss->SetMaterial({{1, 1, 1}}); + +auto fitter = std::make_unique(charge_p, mass_p, std::move(eloss)); +fitter->SetBField({0, 0, 2.85}); +fitter->SetUKFParameters(1e-3, 2, 0); +fitter->SetMeasurementSigma(1.0); + +auto task = new AtFitterTask(std::move(fitter)); +task->SetPersistence(kTRUE); +run->AddTask(task); +``` + +--- + +## MCFitter Assessment + +`AtMCFitter` uses stochastic multi-iteration Monte Carlo optimization and produces `AtMCResult` arrays — not per-track `AtFittedTrack` objects. Its `Exec(AtPatternEvent&)` pattern is fundamentally different from `GetFittedTrack(AtTrack*)`. It does **not** combine well with the new `EventFit::AtFitter` base and should remain separate (`AtMCFitterTask`). + +--- + +## Verification + +1. **Build**: `cmake --build build -j10` — must compile without errors or clang-tidy warnings. +2. **Format**: `clang-format-17 -i` on all new files (or `scripts/formatAll.sh`). +3. **Unit tests**: `cd build && ctest -V -R OpenKFTests` — all new gtest cases must pass. +4. **Existing macro**: Run `macro/tests/UKF/UKFSingleTrack.C` to confirm `TrackFitterUKF` still works standalone. +5. **New macro**: Run `macro/tests/UKF/UKFTask.C` end-to-end on a simulation file; inspect the `AtTrackingEvent` branch in the output ROOT file — verify `AtFittedTrack` contains non-trivial kinematics and smoothed positions. +6. **Check merge**: Verify `AtGenfit` still compiles (it uses `AtFitterOld` which is unchanged). diff --git a/.claude/resume-ukf-task-integration.md b/.claude/resume-ukf-task-integration.md new file mode 100644 index 000000000..86911dd2c --- /dev/null +++ b/.claude/resume-ukf-task-integration.md @@ -0,0 +1,38 @@ +# Resume: UKF Task Integration + +## Status +Plan is fully approved and ready to implement. No code has been written yet. + +## Branch +`OpenKF-Claude` (based on `develop`) + +## Plan file +Full implementation plan is at: `.claude/plans/keen-marinating-wind.md` (inside the repo) + +## Summary of what to do (in order) + +1. **Merge `RealAurio/fitterRefactor` into `OpenKF-Claude`** — this PR (#254) refactors `AtFittedTrack` (new nested structs), introduces `EventFit::AtFitter` base class, adds `AtFitMetadata`/`AtFitTrackMetadata`, and deprecates old classes as `*Old`. 36 files changed. + +2. **Add `fSmoothedPositions` to `AtFittedTrack::TrackProperties`** — `std::vector` for the UKF smoothed trajectory. + +3. **Implement `EventFit::AtFitterUKF`** in `AtReconstruction/AtFitter/AtFitterUKF.h/cxx`: + - Inherits `EventFit::AtFitter`, implements `GetFittedTrack(AtTrack*, ...)` + - Momentum seed from circle-fit Brho (same as GenFit): `brho = B*radius/sin(theta)` → `AtKinematics::GetMomFromBrho()` + - Internally owns a `kf::TrackFitterUKF` (already implemented) + - Outputs smoothed vertex state + smoothed positions into `AtFittedTrack` + - Use plain C++ types (`double`/`int`/`bool`), not ROOT typedefs — class is never persisted + +4. **Update CMakeLists + LinkDef**: + - Add `AtFitter/AtFitterUKF.cxx` to `AtReconstruction/CMakeLists.txt` + - `AtReconstructionLinkDef.h`: `EventFit::AtFitterUKF -!;` + - `AtDataLinkDef.h`: verify `AtFittedTrack +;`, `AtFittedTrack::TrackProperties +;`, new metadata classes + +5. **Write gtests** in `AtReconstruction/AtFitter/AtFitterUKFTest.cxx` (registered in existing `OpenKFTests` target). Cluster positions hardcoded (no file access). Key tests: construction, skip-few-clusters, reasonable kinematics, KE within 20% of truth, smoothed positions populated, metadata converged, `AtFittedTrack` round-trips. + +6. **Write usage macro** `macro/tests/UKF/UKFTask.C`. + +## Key constraints (from CLAUDE.md) +- Tests must NOT access external files — hardcode data; VMCWORKDIR repo files are OK +- `clang-format-17` all new files +- `cmake --build build -j10` +- `cd build && ctest -V -R OpenKFTests` diff --git a/AGENTS.md b/AGENTS.md new file mode 100644 index 000000000..1b5c6eb08 --- /dev/null +++ b/AGENTS.md @@ -0,0 +1,57 @@ +# AGENTS.md + +ATTPCROOT is a ROOT/FairRoot-based C++ framework for simulation and analysis of Active Target Time Projection Chamber (AT-TPC) detector data. + +## Documentation + +Full developer documentation lives in `docs/development/`: + +| Topic | File | +|-------|------| +| Environment setup | [docs/development/setup.md](../docs/development/setup.md) | +| Build & test commands | [docs/development/building.md](../docs/development/building.md) | +| Module overview | [docs/development/modules.md](../docs/development/modules.md) | +| Simulation pipeline | [docs/development/simulation-pipeline.md](../docs/development/simulation-pipeline.md) | +| Reconstruction pipeline | [docs/development/reconstruction-pipeline.md](../docs/development/reconstruction-pipeline.md) | +| Event generators | [docs/development/generators.md](../docs/development/generators.md) | +| Pulse shape analysis | [docs/development/psa.md](../docs/development/psa.md) | +| Energy loss | [docs/development/energy-loss.md](../docs/development/energy-loss.md) | +| Adding a new module | [docs/development/new-module.md](../docs/development/new-module.md) | +| Code style | [docs/development/code-style.md](../docs/development/code-style.md) | +| Testing | [docs/development/testing.md](../docs/development/testing.md) | + +## Quick Reference: Build & Test + +```bash +source build/config.sh # load environment (do this first) +cmake --build build -j10 # build everything +cd build && ctest -V # run all unit tests +``` + +## Code-Writing Rules + +These rules apply whenever editing or adding C++ code in this repo. + +### LinkDef Streamer Suffixes + +Every class in a `*LinkDef.h` file must use the correct suffix: + +- `ClassName +;` — generates a full I/O streamer. **Only** for classes written to a ROOT file (stored in a `TClonesArray` or `TTree` branch). Examples: `AtHit`, `AtEvent`, `AtRawEvent`, `AtMCPoint`. +- `ClassName -!;` — reflection only, no streamer. Use for **everything else**: tasks, algorithms, models, samplers. Examples: `AtPSAMax`, `AtELossModel`, `AtFitterTask`. + +Default to `-!` unless disk persistence is actually required. + +### C++ vs ROOT Typedefs + +- Non-persisted classes (tasks, algorithms, models): use `bool`, `int`, `double`, `std::string`. +- Persisted data classes only: use `Bool_t`, `Int_t`, `Double_t`, etc. + +### Test Isolation + +Unit tests must not access external files or network resources. Hardcode test data inline. + +## Contributing + +- PRs target the `develop` branch; fast-forward only (no merge commits). +- Commit messages: present imperative mood, ≤72 characters. +- All PRs must pass `clang-format`, `clang-tidy`, and unit tests. diff --git a/docs/contributing/code-style.md b/docs/contributing/code-style.md new file mode 100644 index 000000000..33438df3a --- /dev/null +++ b/docs/contributing/code-style.md @@ -0,0 +1,35 @@ +# Code Style + +This page covers formatting and static-analysis tools only. Contributor-critical ROOT and persistence rules live in [guide.md](guide.md). + +## Formatting + +The project uses `clang-format-17` with the repo-root `.clang-format`. + +Format one file: + +```bash +clang-format-17 -i path/to/file.cxx +``` + +Format changed files: + +```bash +scripts/formatAll.sh +``` + +## Static Analysis + +The repo also carries `.clang-tidy` configuration and CI checks for: + +- `clang-format` +- `clang-tidy` +- `iwyu` + +Run clang-tidy manually on a file: + +```bash +clang-tidy path/to/file.cxx -- -I... +``` + +For framework-specific contributor rules, use [guide.md](guide.md). diff --git a/docs/contributing/guide.md b/docs/contributing/guide.md new file mode 100644 index 000000000..73ef7a25a --- /dev/null +++ b/docs/contributing/guide.md @@ -0,0 +1,94 @@ +# Contributor Guide + +## Persisted vs Non-Persisted Code + +Use plain C++ types in non-persisted classes such as tasks, algorithms, models, samplers, and helpers: + +```cpp +bool flag; +int count; +double value; +std::string name; +``` + +Use ROOT typedefs only in classes written to ROOT files: + +```cpp +Bool_t fFlag; +Int_t fCount; +Double_t fValue; +``` + +## ROOT Dictionaries + +Every class with `ClassDef` must appear in the appropriate `*LinkDef.h`. + +Preferred suffix rules for new or edited code: + +- `ClassName +;` + for classes written to disk +- `ClassName -!;` + for reflection-only classes that should not get streamers + +Default to `-!` unless the class is actually persisted. + +If you add a nested type inside a persisted class and ROOT needs to stream it, register that nested type explicitly as well. +If a persisted class owns streamed pointees through pointers or smart pointers, the pointee types still need the required dictionary coverage. + +Current LinkDefs in the tree still contain legacy `+` entries for some non-persisted classes. Treat the rule above as the preferred direction for new work, not as permission for mechanical repo-wide rewrites. + +## Persisted Data Model Conventions + +Do not assume older ROOT folklore applies everywhere in this tree. Persisted classes here already use STL containers and `std::unique_ptr` where ROOT dictionary support is in place. + +Examples in the current branch: + +- `AtEvent` stores hits in `std::vector>` +- `AtTrack` stores hits and pattern state with smart pointers +- `AtTrackingEvent` stores fitted tracks in `std::vector>` + +For new code, follow the local data-model pattern of the owning module. Do not "simplify" these classes back to raw pointers or `TClonesArray` just because older ROOT examples elsewhere do that. + +## Include and Header Norms + +- Public headers are copied into `build/include`, so code in this tree normally includes project headers by basename such as `"AtEvent.h"`, not module-qualified paths. +- Prefer the include style already used in the surrounding module. Avoid mechanical rewrites to `AtData/AtEvent.h`-style paths. +- Public headers commonly use forward declarations for ROOT and FairRoot types. Follow the local header pattern instead of eagerly adding heavyweight includes. + +## When You Add Code + +For a new persisted class, usually update all of: + +- header/source in the owning module +- module `CMakeLists.txt` +- module `*LinkDef.h` +- tests +- docs if branch flow or the runtime data model changes + +For a new non-persisted task or algorithm, usually update all of: + +- header/source in the owning module +- module `CMakeLists.txt` +- module `*LinkDef.h` with `-!` +- tests +- [branch-io-contracts.md](../reference/branch-io-contracts.md) if it reads or writes FairRoot branches + +## CMake and Test Registration + +- `generate_target_and_root_library(...)` + wires libraries and dictionaries together +- `attpcroot_generate_tests(...)` + wires test binaries into CTest + +If you update code but miss one of those registrations, the branch will often build partially and fail later. + +For test writing rules and registration patterns see [testing.md](../tooling/testing.md). + +## ROOT / FairRoot Gotchas + +- Objects read from `FairRootManager` or pulled out of framework-managed `TClonesArray` containers are framework-owned. +- Do not manually `delete` branch objects returned by the framework. +- Do not wrap FairRoot-owned branch objects in smart pointers. +- Macros still use raw `new` for tasks and generators because ownership is transferred through framework registration such as `AddTask(...)` and generator setup. +- Inside tasks and helper classes, prefer the ownership pattern already used locally. In current code that often means `std::unique_ptr` for owned algorithms or models. +- If a ROOT macro cannot resolve a symbol, first suspect the library, dictionary, or environment setup before adding project header includes. diff --git a/docs/contributing/new-module.md b/docs/contributing/new-module.md new file mode 100644 index 000000000..c5c263ad4 --- /dev/null +++ b/docs/contributing/new-module.md @@ -0,0 +1,71 @@ +# Adding a New Module Library + +Each ATTPCROOT module is a separate CMake shared library with an associated ROOT dictionary. + +## CMakeLists.txt Pattern + +```cmake +set(LIBRARY_NAME MyModule) +set(SRCS + MyClass.cxx + AnotherClass.cxx +) +set(DEPENDENCIES + ATTPCROOT::AtData + ROOT::Core + FairRoot::Base +) + +# Optional: unit tests +set(TEST_SRCS MyClassTest.cxx) +attpcroot_generate_tests(${LIBRARY_NAME}Tests SRCS ${TEST_SRCS} DEPS ${LIBRARY_NAME}) + +generate_target_and_root_library( + ${LIBRARY_NAME} + LINKDEF ${LIBRARY_NAME}LinkDef.h + SRCS ${SRCS} + DEPS_PUBLIC ${DEPENDENCIES} +) +``` + +## LinkDef File + +Every class with a `ClassDef` macro in its header must appear in `MyModuleLinkDef.h`. ROOT uses this file to generate the dictionary. + +```cpp +#ifdef __CINT__ +#pragma link off all globals; +#pragma link off all classes; +#pragma link off all functions; + +#pragma link C++ class MyDataClass+; // preferred for persisted classes +#pragma link C++ class MyTask-!; // preferred for reflection-only code +#pragma link C++ class MyAlgorithm-!; + +#endif +``` + +For the preferred streamer suffix rules, persisted-vs-non-persisted type rules, and notes on legacy LinkDef entries already in the tree, use [guide.md](guide.md). + +## Common FairTask Pattern + +Most framework tasks follow: + +- `Init()` + load/register branches and initialize owned algorithms +- `SetParContainers()` + load runtime-db containers when needed +- `Exec()` + clear outputs, read the input `TClonesArray` and current event object, and write the next branch container + +If a task reads or writes documented branches, update [branch-io-contracts.md](../reference/branch-io-contracts.md). + +## Registering with the Build + +Add the new subdirectory to the top-level `CMakeLists.txt`: + +```cmake +add_subdirectory(MyModule) +``` + +And add `ATTPCROOT::MyModule` to the `DEPENDENCIES` of any module that needs it. diff --git a/docs/development/fitting-status.md b/docs/development/fitting-status.md new file mode 100644 index 000000000..bdcdbd78e --- /dev/null +++ b/docs/development/fitting-status.md @@ -0,0 +1,19 @@ +# Fitting Status + +Fitting is unstable on this branch. This page records only current status. + +## Current State + +The current branch contains multiple fitting approaches, including: + +- GenFit-based fitting through `AtFitterTask` and `AtFITTER::AtGenfit` +- Monte-Carlo-based fitting through `AtMCFitterTask` +- Bragg-curve-related fitting logic in other parts of the codebase + +No stable long-term fitting architecture is documented here. + +## Guidance + +- Treat the current branch code as the source of truth. +- Do not assume the current fitters represent the intended long-term design. +- Inspect source directly before editing fitter code. diff --git a/docs/index.md b/docs/index.md new file mode 100644 index 000000000..2a6de11ac --- /dev/null +++ b/docs/index.md @@ -0,0 +1,34 @@ +# Development Documentation + +Agent-facing docs for the current branch state of ATTPCROOT. + +## Tasks + +- install or first-time configure: [tooling/installation.md](tooling/installation.md) +- set up a session / build / run tests: [tooling/daily-use.md](tooling/daily-use.md) +- write and register tests: [tooling/testing.md](tooling/testing.md) +- add code safely: [contributing/guide.md](contributing/guide.md), [contributing/new-module.md](contributing/new-module.md) +- check formatting rules: [contributing/code-style.md](contributing/code-style.md) +- trace branch flow: [reference/branch-io-contracts.md](reference/branch-io-contracts.md), [subsystems/reconstruction-pipeline.md](subsystems/reconstruction-pipeline.md), [subsystems/simulation-pipeline.md](subsystems/simulation-pipeline.md) +- understand runtime objects: [reference/data-model.md](reference/data-model.md) +- find subsystem ownership: [reference/modules.md](reference/modules.md) +- find an example macro: [reference/macro-cookbook.md](reference/macro-cookbook.md) + +## Architecture + +- module map: [reference/modules.md](reference/modules.md) +- data model: [reference/data-model.md](reference/data-model.md) +- simulation flow: [subsystems/simulation-pipeline.md](subsystems/simulation-pipeline.md) +- reconstruction flow: [subsystems/reconstruction-pipeline.md](subsystems/reconstruction-pipeline.md) +- generators: [subsystems/generators.md](subsystems/generators.md) +- PSA: [subsystems/psa.md](subsystems/psa.md) +- energy loss: [subsystems/energy-loss.md](subsystems/energy-loss.md) + +## Branch-Specific Notes + +- fitting status: [fitting-status.md](fitting-status.md) + +## Scope Notes + +- These docs describe the current branch state. +- They are framework-centric and keep experiment-specific detail light. diff --git a/docs/reference/branch-io-contracts.md b/docs/reference/branch-io-contracts.md new file mode 100644 index 000000000..a1d69c09c --- /dev/null +++ b/docs/reference/branch-io-contracts.md @@ -0,0 +1,33 @@ +# Branch and IO Contracts + +This page records the current branch names and persistence defaults used by major tasks. These defaults come from the task constructors and `Init()` implementations in the current branch. + +In most cases the branch object returned by `FairRootManager::GetObject(...)` is a `TClonesArray`, not a bare event object. The task then reads slot `0` as the current `AtRawEvent`, `AtEvent`, `AtPatternEvent`, or `AtTrackingEvent`. + +## Reconstruction Tasks + +| Task | Input Branches | Output Branches | Persistence Default | Notes | +|------|----------------|-----------------|---------------------|-------| +| `AtUnpackTask` | external file via unpacker | `AtRawEvent` | `true` | Registers a `TClonesArray` holding `AtRawEvent` objects | +| `AtFilterTask` | `AtRawEvent` | `AtRawEventFiltered` | `false` | Reads/writes `TClonesArray`-style branches; can also filter aux/FPN pads | +| `AtPSAtask` | `AtRawEvent` | `AtEventH` | `false` | Reads `TClonesArray` input, writes `TClonesArray("AtEvent", 1)` output; optionally looks for `AtTpcPoint` | +| `AtDataCleaningTask` | `AtEventH` | `AtEventCleaned` | `true` | Reads the input event from a `TClonesArray` and writes a cleaned `AtEvent` copy | +| `AtPRAtask` | `AtEventH` | `AtPatternEvent` | `false` | PRA-based pattern recognition path using `TClonesArray` event containers | +| `AtSampleConsensusTask` | `AtEventH` | `AtPatternEvent` | `false` | Sample-consensus path using the same container pattern | +| `AtFitterTask` | `AtPatternEvent` | `AtTrackingEvent` | `false` | Output is a `TClonesArray("AtTrackingEvent", 1)`; verify input-branch behavior in source before relying on setters | +| `AtMCFitterTask` | `AtPatternEvent` | `AtMCResult`, `SimEvent`, `SimRawEvent` | `true`, `false`, `false` | Save toggles exist for each output branch | + +## Simulation and Digitization Tasks + +| Task | Input Branches | Output Branches | Persistence Default | Notes | +|------|----------------|-----------------|---------------------|-------| +| `AtClusterizeTask` | `AtTpcPoint` | `AtSimulatedPoint` | `false` | Reads the detector-specific `AtTpcPoint` `TClonesArray`, not a generic `AtMCPoint` branch | +| `AtPulseTask` | `AtSimulatedPoint`, optional `AtTpcPoint` | `AtRawEvent`, optional re-registered `AtTpcPoint` | `true`, `false` | Reads branch containers, not bare objects; `SetSaveMCInfo()` enables MC mapping behavior | + +## Reading This Page + +- These are default names, not hard-coded universal laws. +- Many tasks expose `SetInputBranch(...)`, `SetOutputBranch(...)`, or related setters. +- Do not assume every setter is honored consistently; check the task source before automating branch-name rewrites. +- If you add a task that reads/writes FairRoot branches, update this page. +- For the runtime objects behind these branches, see [data-model.md](data-model.md). diff --git a/docs/reference/data-model.md b/docs/reference/data-model.md new file mode 100644 index 000000000..a40a1fb6a --- /dev/null +++ b/docs/reference/data-model.md @@ -0,0 +1,46 @@ +# Data Model + +Core runtime objects passed between simulation, unpacking, reconstruction, and fitting. + +Most FairRoot branches in this tree are not bare `AtEvent*` or `AtRawEvent*` objects. They are `TClonesArray` branch containers, usually holding one event object at slot `0`. Agents should treat the branch name and the contained event class as related but distinct concepts. + +## Reconstruction Objects + +| Object | Typical branch container | Produced by | Consumed by | Persisted | Key downstream fields | +|--------|--------------------------|-------------|-------------|-----------|-----------------------| +| `AtRawEvent` | `TClonesArray` branch such as `AtRawEvent` / `AtRawEventFiltered` | `AtUnpackTask`, `AtPulseTask` | filters, PSA | yes | pad traces, aux/FPN pads, good flag, optional MC map | +| `AtEvent` | `TClonesArray` branch such as `AtEventH` / `AtEventCleaned` | `AtPSAtask` | cleaning, pattern recognition | yes | hit list, event charge, mesh signal | +| `AtPatternEvent` | `TClonesArray` branch `AtPatternEvent` | `AtPRAtask`, `AtSampleConsensusTask` | fitting | yes | candidate tracks, noise hits | +| `AtTrack` | contained inside `AtPatternEvent`, not its own top-level branch | pattern recognition code | fitters | yes | track ID, hit/cluster-hit collections, geometric estimates, Bragg-curve values | +| `AtTrackingEvent` | `TClonesArray` branch `AtTrackingEvent` | `AtFitterTask` | downstream analysis | yes | fitted tracks, optional copied track array, event vertex fields | +| `AtFittedTrack` | contained inside `AtTrackingEvent`, not its own top-level branch | fitter implementations | downstream analysis | yes | current branch layout only; treat as unstable | + +## Simulation Objects + +| Object | Typical branch container | Produced by | Consumed by | Persisted | Notes | +|--------|--------------------------|-------------|-------------|-----------|-------| +| `AtTpcPoint` | `TClonesArray` branch `AtTpcPoint` | Geant/VMC transport through the AT-TPC detector | digitization tasks such as `AtClusterizeTask` and `AtPulseTask` | yes | detector-specific point branch seen by downstream digitization code | +| `AtMCPoint` | module-specific MC point branches | Geant/VMC transport | analysis, truth matching, detector-specific point types | yes | generic MC-point base family; downstream digitization usually works with detector-specific derived point classes | +| `AtMCTrack` | `TClonesArray` branch `MCTrack` | simulation stack | analysis, truth matching | yes | simulated particle tracks | +| `AtVertexPropagator` | no FairRoot branch | generators/runtime simulation logic | generators and downstream simulation code | no | singleton shared state | + +## Flow Summary + +```text +TClonesArray(AtRawEvent branch) -> AtRawEvent +TClonesArray(AtEventH branch) -> AtEvent +TClonesArray(AtPatternEvent) -> AtPatternEvent +TClonesArray(AtTrackingEvent) -> AtTrackingEvent + | | + | -> contains AtTrack + -> contains hits +``` + +```text +Generator / AtVertexPropagator + -> MCTrack / AtTpcPoint branches + -> AtSimulatedPoint branch + -> AtRawEvent branch +``` + +For the tasks that produce and consume these objects, see [branch-io-contracts.md](branch-io-contracts.md). For pipeline order, see [simulation-pipeline.md](../subsystems/simulation-pipeline.md) and [reconstruction-pipeline.md](../subsystems/reconstruction-pipeline.md). diff --git a/docs/reference/macro-cookbook.md b/docs/reference/macro-cookbook.md new file mode 100644 index 000000000..0ff481383 --- /dev/null +++ b/docs/reference/macro-cookbook.md @@ -0,0 +1,26 @@ +# Macro Cookbook + +Example macro starting points referenced by these docs. This page is intentionally sparse. + +## Notes + +- `macro/examples/` contains adaptation starting points. +- `macro/tests/` contains repeatable workflow and regression-style examples. +- Many macros are experiment-specific or user-specific and are intentionally omitted. +- Files under `macro/` are ROOT / Cling scripts, not normal compiled translation units. +- ROOT/core headers such as `TClonesArray.h` or `TTreeReader.h`, and STL headers a macro genuinely uses, are normal. +- Do not treat a macro like standalone C++ or try to fix missing symbols by adding project headers first. + +## Macro Do / Don't + +- Do: copy an existing pattern from `macro/examples/` or `macro/tests/`. +- Don't: start by adding `#include "At*.h"` to compensate for missing dictionaries or libraries, or rewrite a macro as if it were a compiled source file. + +## Starting Points + +| Path | Use | Caveat | +|------|-----|--------| +| `macro/examples/run_sim.C` | minimal simulation structure | simple example, not a full workflow template | +| `macro/examples/rundigi_sim.C` | simulation plus downstream task ordering | parts of the task usage are older and should be checked against current APIs | +| `macro/tests/AT-TPC/run_unpack_attpc.C` | unpacking plus PSA/cleaning flow | experiment-facing example, not a framework contract | +| `macro/tests/AT-TPC/run_sim_attpc.C` | fuller simulation test flow | oriented toward test workflow rather than minimal onboarding | diff --git a/docs/reference/modules.md b/docs/reference/modules.md new file mode 100644 index 000000000..4182e3228 --- /dev/null +++ b/docs/reference/modules.md @@ -0,0 +1,28 @@ +# Module Overview + +ATTPCROOT is organized as CMake library targets layered around data, simulation, unpacking, reconstruction, and analysis. + +## Modules + +| Module | Role | Main subareas | Common edit reasons | Also inspect | +|--------|------|---------------|---------------------|--------------| +| `AtData` | persisted event and track containers | event, track, pattern classes | add/modify data objects | LinkDef, producing tasks, [data-model.md](data-model.md) | +| `AtReconstruction` | PSA, filtering, pattern recognition, fitting tasks | `AtPulseAnalyzer/`, `AtPatternRecognition/`, `AtFitter/`, `AtFilter/` | change branch flow or task behavior | [branch-io-contracts.md](branch-io-contracts.md), tests | +| `AtUnpack` | experimental data to `AtRawEvent` | unpackers, `GETDecoder2/` | add file-format/front-end support | maps, macros | +| `AtDigitization` | simulation truth to detector signals | clusterization, pulse generation | change MC-to-signal behavior | simulation pipeline | +| `AtTools` | shared utilities | energy loss, kinematics, hit sampling, cleaning | algorithmic utilities used across modules | dependent modules | +| `AtSimulationData` | MC truth and simulation shared state | MC objects, propagator | change truth object layout or simulation state | simulation pipeline | +| `AtGenerators` | FairRoot generator implementations | reaction and beam generators | change simulation source behavior | `AtVertexPropagator`, simulation pipeline | +| `AtDetectors` | detector geometry/sensitive detectors | detector-specific subdirs | detector implementation changes | geometry assets | +| `AtAnalysis` | higher-level analysis code | experiment-facing analysis | analysis task changes | data model | +| `AtEventDisplay` | event visualization | tabs/display code | visualization changes | data model | + +## Workflow Assets + +| Directory | Role | +|-----------|------| +| `macro/` | example and test macros | +| `geometry/` | detector geometry assets | +| `parameters/` | runtime parameter files | +| `resources/` | shipped data assets such as tables and cross sections | +| `scripts/` | support scripts and mapping-related assets | diff --git a/docs/subsystems/energy-loss.md b/docs/subsystems/energy-loss.md new file mode 100644 index 000000000..1f1b419af --- /dev/null +++ b/docs/subsystems/energy-loss.md @@ -0,0 +1,59 @@ +# Energy Loss + +ATTPCROOT provides several energy-loss utilities in `AtTools/`. The modern model interface is `AtTools::AtELossModel`; `AtTools::AtELossManager` is a separate legacy lookup-table helper, not the owner or facade for all models. + +## Main Types + +| Type | Role | +|------|------| +| `AtELossModel` | abstract interface for stopping power, range, energy loss, and residual energy | +| `AtELossCATIMA` | CATIMA-backed implementation | +| `AtELossTable` | table-backed implementation, typically from SRIM-style data | +| `AtELossBetheBloch` | standalone Bethe-Bloch helper for simpler analytic use cases | +| `AtELossManager` | older lookup-table utility; do not treat it as the generic `AtELossModel` entry point | + +## `AtELossModel` Interface + +All `AtELossModel` implementations expose the same core API: + +- `GetdEdx(energy)` +- `GetRange(energyIni, energyFin = 0)` +- `GetEnergyLoss(energyIni, distance)` +- `GetEnergy(energyIni, distance)` + +Units in this layer are MeV and mm unless a class-specific note says otherwise. + +## CATIMA + +`AtELossCATIMA` is the main model-backed implementation. It wraps the CATIMA library and requires both material configuration and projectile configuration. + +Typical setup: + +```cpp +using MaterialComp = std::tuple; + +std::vector material = { + {1, 1, 2}, // example stoichiometry entries + {12, 6, 1}, +}; + +AtTools::AtELossCATIMA model(density_g_cm3, material); +model.SetProjectile(projectileA, projectileZ, projectileMassAmu); + +double dE = model.GetEnergyLoss(energyMeV, distanceMm); +double range = model.GetRange(energyMeV); +double residual = model.GetEnergy(energyMeV, distanceMm); +``` + +CATIMA is fetched by CMake if it is not already available locally. It is the only backend in this tree that currently exposes straggling calculations through the `AtELossModel` interface. + +## Table and Analytic Helpers + +- `AtELossTable` follows the same `AtELossModel` interface, but uses precomputed stopping/range tables. +- `AtELossBetheBloch` is a lighter analytic helper and is not a replacement for the full CATIMA-backed workflow when you need straggling or richer material handling. + +## Legacy Lookup-Table Path + +`AtELossManager` predates the `AtELossModel` hierarchy. It still exists in the tree, but it is a separate lookup-table class with its own API such as `GetEnergyLoss(...)`, `GetFinalEnergy(...)`, and `GetDistance(...)`. + +Use it when you specifically need that legacy behavior; do not document or refactor it as if it were the generic model manager for the rest of `AtTools/`. diff --git a/docs/subsystems/generators.md b/docs/subsystems/generators.md new file mode 100644 index 000000000..289d3eb3a --- /dev/null +++ b/docs/subsystems/generators.md @@ -0,0 +1,36 @@ +# Event Generators + +ATTPCROOT uses FairRoot's generator framework. All reaction generators subclass `AtReactionGenerator`, which itself subclasses `FairGenerator`. + +## AtReactionGenerator + +`AtReactionGenerator` is the abstract base for all physics reaction generators in `AtGenerators/`. + +**What to implement:** Subclasses implement `GenerateReaction(FairPrimaryGenerator*)`. This method is responsible for computing the kinematics of the reaction and adding primaries to the stack. + +**What not to touch:** `ReadEvent()` is declared `final`. It handles the beam/reaction event alternation via `AtVertexPropagator` and calls `GenerateReaction()` internally. Do not override it. + +## AtVertexPropagator + +`AtVertexPropagator` is a singleton that passes vertex position, beam momentum, and kinematics between generators and downstream tasks (e.g., digitization tasks that need beam parameters). + +```cpp +AtVertexPropagator::Instance()->GetVertex(); // read vertex +AtVertexPropagator::Instance()->EndEvent(); // called at end of each event to alternate beam/reaction +``` + +**In unit tests:** Call `AtVertexPropagator::Instance()->ResetForTesting()` in `SetUp()` to start from a clean state between tests. + +## Sequential Decays + +Generators can be chained for multi-step decays. Enable on the upstream generator: + +```cpp +generator->SetSequentialDecay(true); +``` + +The second generator in the chain reads the vertex and momentum from `AtVertexPropagator` rather than sampling a fresh beam. + +## Beam/Reaction Alternation + +By default, `ReadEvent()` alternates between inserting a beam event (no reaction) and a full reaction event. This models the AT-TPC's operation where beam tracks are recorded alongside reaction tracks for calibration. The alternation is controlled via `AtVertexPropagator::EndEvent()`. diff --git a/docs/subsystems/psa.md b/docs/subsystems/psa.md new file mode 100644 index 000000000..2cfae634a --- /dev/null +++ b/docs/subsystems/psa.md @@ -0,0 +1,40 @@ +# Pulse Shape Analysis (PSA) + +PSA extracts hit positions and charge information from raw pad traces (`AtRawEvent`) to produce an `AtEvent` containing `AtHit` objects. + +## Class Hierarchy + +`AtPSA` is the abstract base class in `AtReconstruction/AtPulseAnalyzer/`. All concrete PSA implementations inherit from it. + +| Class | Algorithm | +|-------|-----------| +| `AtPSAMax` | Peak-finding: hit at the time bucket of maximum signal | +| `AtPSAFull` | Integrates the full trace; produces one hit per pad | +| `AtPSADeconv` | Deconvolution-based hit finding; can find multiple hits per pad | + +## Running PSA + +PSA is applied by `AtPSAtask`, which wraps any `AtPSA` subclass: + +```cpp +auto psa = std::make_unique(); +auto psaTask = new AtPSAtask(std::move(psa)); +fRun->AddTask(psaTask); +``` + +## Chaining PSA Methods + +`AtPSAComposite` allows running multiple PSA methods in sequence on the same event, for example to handle different pad types with different algorithms: + +```cpp +auto composite = std::make_unique(); +composite->AddPSA(std::make_unique()); +composite->AddPSA(std::make_unique()); +``` + +## Output + +PSA produces an `AtEvent` stored as a branch in the output TTree. Each `AtHit` in the event carries: +- 3D position (x, y, time-bucket → z via drift velocity) +- Integrated charge +- Pad number and other metadata diff --git a/docs/subsystems/reconstruction-pipeline.md b/docs/subsystems/reconstruction-pipeline.md new file mode 100644 index 000000000..4845df473 --- /dev/null +++ b/docs/subsystems/reconstruction-pipeline.md @@ -0,0 +1,62 @@ +# Reconstruction Pipeline + +The reconstruction pipeline converts raw pad traces into hit clouds, candidate tracks, and fitted track results. Each stage is usually a `FairTask` that reads and writes branches through `FairRootManager`. + +In the current tree, those branches are usually `TClonesArray` containers. Tasks typically fetch the branch container from `FairRootManager`, then read or construct the event object at slot `0`. + +## Shared Early Stages + +```text +AtUnpackTask + │ reads raw GET/GRAW, HDF5, or ROOT data + └─ output branch: AtRawEvent -> TClonesArray[AtRawEvent] + ▼ +AtFilterTask [optional] + │ filters raw traces + └─ output branch: AtRawEventFiltered -> TClonesArray[AtRawEvent] + ▼ +AtPSAtask + │ converts traces into hits + └─ output branch: AtEventH -> TClonesArray[AtEvent] + ▼ +AtDataCleaningTask [optional] + │ removes or adjusts hits before pattern recognition + └─ output branch: AtEventCleaned -> TClonesArray[AtEvent] +``` + +The exact branch names and defaults are documented in [branch-io-contracts.md](../reference/branch-io-contracts.md). + +## Pattern Recognition + +Two supported paths currently produce `AtPatternEvent`: + +```text +AtPRAtask + │ PRA / triplclust-style pattern recognition + └─ output branch: AtPatternEvent -> TClonesArray[AtPatternEvent] +``` + +```text +AtSampleConsensusTask + │ sample-consensus pattern recognition + └─ output branch: AtPatternEvent -> TClonesArray[AtPatternEvent] +``` + +These are both supported. The docs do not declare one of them to be the universal default. + +## Fitting + +After pattern recognition, fitting produces the `AtTrackingEvent` branch container. + +Current fit-related tasks include: + +- `AtFitterTask` +- `AtMCFitterTask` + +Fitting is unstable, so this page intentionally stays shallow. See [fitting-status.md](../development/fitting-status.md). + +## Running Reconstruction + +A typical macro creates a `FairRunAna`, adds tasks in order, and calls `Run(...)`. + +For runtime objects, see [data-model.md](../reference/data-model.md). For example macros, see [macro-cookbook.md](../reference/macro-cookbook.md). diff --git a/docs/subsystems/simulation-pipeline.md b/docs/subsystems/simulation-pipeline.md new file mode 100644 index 000000000..cccb33414 --- /dev/null +++ b/docs/subsystems/simulation-pipeline.md @@ -0,0 +1,54 @@ +# Simulation Pipeline + +The simulation pipeline converts a generated reaction into MC truth, simulated detector responses, and finally `AtRawEvent` output that can be passed into reconstruction. + +Like reconstruction, the task-facing branch objects here are usually `TClonesArray` containers. The detector-point branch used by digitization is `AtTpcPoint`, not a generic `AtMCPoint` branch. + +## Steps + +``` +FairPrimaryGenerator + └─ AtReactionGenerator subclass + │ generates primary particles and writes to particle stack + ▼ +Geant4 / VMC transport + │ produces detector-point branches such as AtTpcPoint + ▼ +AtClusterizeTask + │ reads TClonesArray(AtTpcPoint) and converts deposits -> ionization electron clusters + ▼ +AtPulseTask + │ drifts electrons and produces simulated pad traces + └─ output branch: AtRawEvent -> TClonesArray[AtRawEvent] +``` + +## Key Runtime Objects + +- `AtVertexPropagator` + shared simulation-side state between generators and downstream logic +- `AtMCTrack` + simulated particle tracks +- `AtTpcPoint` + detector-specific transport hits consumed by current digitization tasks +- `AtMCPoint` + generic base-family concept for detector points, but usually not the concrete branch a digitization task reads +- `AtRawEvent` + final simulated raw traces used by reconstruction + +See [data-model.md](../reference/data-model.md) for the object-level view and [branch-io-contracts.md](../reference/branch-io-contracts.md) for task branch names. + +## Running a Simulation + +Simulation macros live mostly in `macro/Simulation/` and `macro/examples/`. A typical macro: + +1. Creates a `FairRunSim` +2. Attaches the detector geometry +3. Configures a `FairPrimaryGenerator` with one or more `AtReactionGenerator` subclasses +4. Adds `AtClusterizeTask` and `AtPulseTask` +5. Calls `fRun->Run(nEvents)` + +The output ROOT file contains MC truth alongside detector-specific point branches and the simulated `AtRawEvent`, so it can be fed directly into the reconstruction pipeline. + +Geometry comes from `geometry/`; runtime parameters typically come from `parameters/`. + +See [generators.md](generators.md) for generator behavior and [macro-cookbook.md](../reference/macro-cookbook.md) for example macros. diff --git a/docs/tooling/daily-use.md b/docs/tooling/daily-use.md new file mode 100644 index 000000000..8746256b8 --- /dev/null +++ b/docs/tooling/daily-use.md @@ -0,0 +1,61 @@ +# Daily Use + +## Before Every Session + +Source the generated environment script from the repo root: + +```bash +source build/config.sh +``` + +This sets the following environment variables: + +| Variable | Purpose | +|----------|---------| +| `ROOTSYS` | ROOT installation prefix | +| `VMCWORKDIR` | Repo root (used by macros to find geometry, parameters, etc.) | +| `ROOT_INCLUDE_PATH` | Headers exposed to the ROOT interpreter | +| `LD_LIBRARY_PATH` | Runtime library search path | +| Geant4 data paths | `G4...DATA` variables needed by Geant4 | + +VSCode terminals auto-source this via the workspace settings. + +If `build/config.sh` does not exist, the project has not been configured yet — see [installation.md](installation.md). + +## Building + +```bash +cmake --build build -j10 +``` + +Adjust `-j` to your CPU count. To rebuild a single module: + +```bash +cmake --build build --target AtSimulationData -j10 +``` + +## Running Tests + +Unit tests are built by default (`BUILD_TESTS=ON`). Binaries land in `build/tests/`. + +Run all tests: + +```bash +cd build && ctest -V +``` + +Run a specific test binary directly: + +```bash +./build/tests/AtSimulationDataTests +./build/tests/AtGeneratorsTests +./build/tests/AtToolsTests +``` + +See [testing.md](testing.md) for how to write and register new tests. + +## Verifying the Environment + +```bash +root-config --version # prints ROOT version +``` diff --git a/docs/tooling/installation.md b/docs/tooling/installation.md new file mode 100644 index 000000000..94c6091f4 --- /dev/null +++ b/docs/tooling/installation.md @@ -0,0 +1,50 @@ +# Installation (First-Time Setup) + +This page covers the one-time steps to configure and build ATTPCROOT on a new machine. After completing these steps, day-to-day work uses only `source build/config.sh` and `cmake --build build` — the environment variables below are never needed again. + +## Prerequisites + +FairSoft and FairRoot must already be installed on the machine. FairSoft provides ROOT, Geant4, and VMC. FairRoot provides the task/run framework. + +## Step 1 — Set Configure-Time Variables + +These two variables tell CMake where FairRoot and FairSoft are installed. They only need to be set in the shell session where you run `cmake` for the first time: + +```bash +export FAIRROOTPATH=/path/to/your/fairroot/install +export SIMPATH=/path/to/your/fairsoft/install +``` + +These are machine-specific. Ask your sysadmin or check the FairRoot/FairSoft install documentation for the correct paths on your system. + +## Step 2 — Configure + +Run CMake out-of-source from the repo root: + +```bash +cmake -S . -B build +``` + +If HDF5 is not bundled with your FairSoft installation, add its prefix: + +```bash +cmake -S . -B build -DCMAKE_PREFIX_PATH=/path/to/hdf5 +``` + +A successful configure generates `build/config.sh`. This file captures the full runtime environment and is sourced before any subsequent build or analysis session. + +## Step 3 — Build + +```bash +cmake --build build -j10 +``` + +Adjust `-j` to your CPU count. + +## After Installation + +From this point on: + +- **Before any session and to rebuild:** see [daily-use.md](daily-use.md) + +The `FAIRROOTPATH` and `SIMPATH` variables from Step 1 are baked into `build/config.sh` and do not need to be set again. diff --git a/docs/tooling/testing.md b/docs/tooling/testing.md new file mode 100644 index 000000000..be5a8a78e --- /dev/null +++ b/docs/tooling/testing.md @@ -0,0 +1,59 @@ +# Testing + +ATTPCROOT has two categories of tests: unit tests (compiled, fast) and integration tests (interpreted ROOT macros). + +## Unit Tests + +Unit tests use [Google Test](https://github.com/google/googletest). They live alongside the source files they test — for example, `AtSimulationData/AtVertexPropagatorTest.cxx` tests `AtVertexPropagator`. + +### Writing a Test + +Create a `*Test.cxx` file in the module directory and register it in the module's `CMakeLists.txt`: + +```cmake +set(TEST_SRCS MyClassTest.cxx) +attpcroot_generate_tests(${LIBRARY_NAME}Tests SRCS ${TEST_SRCS} DEPS ${LIBRARY_NAME}) +``` + +`attpcroot_generate_tests` links against Google Test and registers the binary with CTest. The test binary is placed in `build/tests/`. + +Typical binary names follow the module, for example: + +- `AtDataTests` +- `AtGeneratorsTests` +- `AtToolsTests` +- `OpenKFTests` for Eigen/OpenKF-specific reconstruction tests + +### Rules + +- **Tests must not access external files or network resources.** Hardcode small data values inline or use files that ship with the repo under `$VMCWORKDIR`. Tests that require external state are fragile and break on other machines. +- Tests should be fast. Unit tests that take more than a second are a smell. +- For naming conventions see the [wiki](https://github.com/ATTPC/ATTPCROOTv2/wiki/Coding-conventions#test-names). + +### Running Unit Tests + +```bash +cd build && ctest -V # all tests +./build/tests/AtToolsTests # one binary directly +``` + +## Integration Tests + +More involved tests (those that require data files or test full pipelines) live in `macro/tests/`. They are written as ROOT macros (`.C` files) and run interpreted, not compiled. + +These are listed and run separately from unit tests. Use them for: + +- end-to-end workflow checks +- simulation plus reconstruction flows +- experiment-facing regression scenarios + +See [macro-cookbook.md](macro-cookbook.md) for curated macro starting points. + +## Continuous Integration + +Every PR runs: +1. `clang-format` check +2. `clang-tidy` and `iwyu` static analysis +3. Full unit test build and run (`ctest`) + +Run the tests locally before submitting a PR to catch failures early. From 9eaab3530bcaecbed742054637e3e2b2641f27f4 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Wed, 1 Apr 2026 14:21:27 -0400 Subject: [PATCH 05/10] Fix MeV/cm vs MeV/mm unit error in AtELossCATIMA MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit GetBraggCurve stored result.dEdxi * fDensity directly (MeV/cm) into the output vector, while all other paths (GetdEdx, base-class GetBraggCurve) divide by 10 to convert to MeV/mm. The energy-step calculation in the same loop used dEdx * rangeStepSize / 10, which accidentally cancelled the unit error during tracking, so range and energy results were unaffected — only the returned dEdx y-values were 10× too large, inflating the displayed Bragg peak height relative to BB and SRIM curves. Adds a fourth model per scenario — CATIMA with bare projectile charge (z_eff_type::none) and dEdx corrections disabled (no Barkas, Lindhard, or shell corrections) — to isolate which physics differences drive the disagreement between Bethe-Bloch and full CATIMA: - dEdx/range divergence at low energy → effective charge (Pierce-Blann) - Range straggling floor in BB vs steep fall in CATIMA → Firsov/Lindhard-X straggling correction (hardwired in CATIMA, absent in BB) - Bragg peak shape difference → effective charge suppressing dEdx as v→0 --- AtTools/AtELossCATIMA.cxx | 4 +- macro/tests/ELossComparison/ELossComparison.C | 75 +++++++++++++++++-- 2 files changed, 69 insertions(+), 10 deletions(-) diff --git a/AtTools/AtELossCATIMA.cxx b/AtTools/AtELossCATIMA.cxx index e8b7919d2..2c26e5cd2 100644 --- a/AtTools/AtELossCATIMA.cxx +++ b/AtTools/AtELossCATIMA.cxx @@ -143,8 +143,8 @@ AtELossCATIMA::GetBraggCurve(double energy, double rangeStepSize, double totalFr while (remainingEnergy / energy > totalFractionELoss) { catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); - double dEdx = result.dEdxi * fDensity; - braggCurve.push_back(std::make_pair(dEdx, range)); + double dEdx = result.dEdxi * fDensity; // MeV/cm + braggCurve.push_back(std::make_pair(dEdx / 10.0, range)); // store in MeV/mm double DE = dEdx * rangeStepSize / 10.; if (DE > remainingEnergy) diff --git a/macro/tests/ELossComparison/ELossComparison.C b/macro/tests/ELossComparison/ELossComparison.C index ddc91a742..43b19be63 100644 --- a/macro/tests/ELossComparison/ELossComparison.C +++ b/macro/tests/ELossComparison/ELossComparison.C @@ -26,6 +26,29 @@ * root -l -b -q ELossComparison.C * * Output: ELossComparison.pdf (four pages, one per scenario) + * + * Model notes + * ----------- + * "Bethe-Bloch (Bohr)" — AtELossBetheBloch: pure PDG 2022 Eq. 34.1, bare projectile charge, + * Bohr range-straggling integral (dω²/dE = K·z²·(Z/A)·ρ·mₑ / |dEdx|³), √E extrapolation + * below the BB validity threshold. + * + * "CATIMA (default)" — AtELossCATIMA with catima::default_config: + * - Effective charge: Pierce-Blann (z_eff → 0 for protons at low β — suppresses dEdx and + * straggling proportionally to z_eff²). + * - dEdx: Bethe-Bloch + shell corrections + Barkas term + density effect + Lindhard correction; + * SRIM-85 tables at low energies. + * - Straggling: Bohr × Lindhard-X factor (quantum correction, always active); Firsov formula + * used as upper cap below 30 MeV/u — this is LOWER than pure Bohr and causes the straggling + * to drop steeply toward zero at low energy. + * + * "CATIMA (bare, no corr)" — AtELossCATIMA with custom config: + * - Effective charge: none (bare Z, matching BB). + * - dEdx corrections disabled: no Barkas, no Lindhard, no shell corrections. + * - Straggling: still Bohr × Lindhard-X + Firsov cap (not removable via config in the bundled + * CATIMA version). Low-energy dEdx still from SRIM-85 tables, not √E. + * This is the closest CATIMA can get to a pure BB+Bohr calculation; any remaining gap versus the + * red BB curve is due to the Firsov straggling correction and SRIM-85 vs √E extrapolation. */ // --------------------------------------------------------------------------- @@ -356,12 +379,44 @@ void ELossComparison() auto table_p_He = LoadSimpleSrimTable(Form("%s/proton_He_700torr.txt", resDir.Data()), kRho_He_700); + // ========================================================================= + // "CATIMA bare" config: bare charge + no dEdx corrections — closest analog + // to the pure Bethe-Bloch model available within the CATIMA library. + // Note: the Lindhard-X / Firsov straggling correction is always active in + // the bundled CATIMA source and cannot be disabled via Config. + // ========================================================================= + catima::Config pureCfg; + pureCfg.z_effective = catima::z_eff_type::none; // bare Z, not Pierce-Blann + pureCfg.corrections = + catima::no_barkas | catima::no_lindhard | catima::no_shell_correction; + + auto catimaPure_p_H2 = + new AtTools::AtELossCATIMA(kRho_H2_600, std::vector>{{1, 1, 1}}); + catimaPure_p_H2->SetProjectile(1, 1, kMproton_amu); + catimaPure_p_H2->SetConfig(pureCfg); + + auto catimaPure_d_D2 = + new AtTools::AtELossCATIMA(kRho_D2_600, std::vector>{{2, 1, 1}}); + catimaPure_d_D2->SetProjectile(2, 1, kMdeuteron_amu); + catimaPure_d_D2->SetConfig(pureCfg); + + auto catimaPure_a_D2 = + new AtTools::AtELossCATIMA(kRho_D2_600, std::vector>{{2, 1, 1}}); + catimaPure_a_D2->SetProjectile(4, 2, kMalpha_amu); + catimaPure_a_D2->SetConfig(pureCfg); + + auto catimaPure_p_He = + new AtTools::AtELossCATIMA(kRho_He_700, std::vector>{{4, 2, 1}}); + catimaPure_p_He->SetProjectile(1, 1, kMproton_amu); + catimaPure_p_He->SetConfig(pureCfg); + // ========================================================================= // Color / line-style scheme // ========================================================================= const int kColBB = kRed + 1; const int kColCATIMA = kBlue + 1; const int kColTable = kGreen + 2; + const int kColCATIMApure = kMagenta + 1; // ========================================================================= // Draw and save @@ -369,9 +424,10 @@ void ELossComparison() auto *c1 = new TCanvas("c1", "Proton in H2", 1200, 900); { std::vector ents = { - {bb_p_H2, "Bethe-Bloch", kColBB, 1, true}, - {catima_p_H2, "CATIMA", kColCATIMA, 2, true}, + {bb_p_H2, "Bethe-Bloch (Bohr)", kColBB, 1, true}, + {catima_p_H2, "CATIMA (default)", kColCATIMA, 2, true}, {table_p_H2, "SRIM Table", kColTable, 3, true}, + {catimaPure_p_H2, "CATIMA (bare, no corr)", kColCATIMApure, 7, true}, }; DrawScenario("Proton in H_{2} (600 Torr)", ents, 0.1, 10.0, 5.0, c1); } @@ -379,8 +435,9 @@ void ELossComparison() auto *c2 = new TCanvas("c2", "Deuteron in D2", 1200, 900); { std::vector ents = { - {bb_d_D2, "Bethe-Bloch", kColBB, 1, true}, - {catima_d_D2, "CATIMA", kColCATIMA, 2, true}, + {bb_d_D2, "Bethe-Bloch (Bohr)", kColBB, 1, true}, + {catima_d_D2, "CATIMA (default)", kColCATIMA, 2, true}, + {catimaPure_d_D2, "CATIMA (bare, no corr)", kColCATIMApure, 7, true}, }; if (table_d_D2) ents.push_back({table_d_D2, "SRIM Table", kColTable, 3, false}); @@ -390,8 +447,9 @@ void ELossComparison() auto *c3 = new TCanvas("c3", "Alpha in D2", 1200, 900); { std::vector ents = { - {bb_a_D2, "Bethe-Bloch", kColBB, 1, true}, - {catima_a_D2, "CATIMA", kColCATIMA, 2, true}, + {bb_a_D2, "Bethe-Bloch (Bohr)", kColBB, 1, true}, + {catima_a_D2, "CATIMA (default)", kColCATIMA, 2, true}, + {catimaPure_a_D2, "CATIMA (bare, no corr)", kColCATIMApure, 7, true}, }; if (table_a_D2) ents.push_back({table_a_D2, "SRIM Table", kColTable, 3, false}); @@ -401,8 +459,9 @@ void ELossComparison() auto *c4 = new TCanvas("c4", "Proton in He", 1200, 900); { std::vector ents = { - {bb_p_He, "Bethe-Bloch", kColBB, 1, true}, - {catima_p_He, "CATIMA", kColCATIMA, 2, true}, + {bb_p_He, "Bethe-Bloch (Bohr)", kColBB, 1, true}, + {catima_p_He, "CATIMA (default)", kColCATIMA, 2, true}, + {catimaPure_p_He, "CATIMA (bare, no corr)", kColCATIMApure, 7, true}, }; if (table_p_He) ents.push_back({table_p_He, "SRIM Table", kColTable, 3, false}); From 5dae8c1f8c35080f3d11a60d96c4d8f0a0454fa8 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Wed, 1 Apr 2026 14:37:01 -0400 Subject: [PATCH 06/10] Update checkout github action version Also fixes clang-formatting --- .github/workflows/CI-tests.yml | 2 +- .github/workflows/main.yml | 2 +- AGENTS.md | 27 +++++++++++++++------------ AtTools/AtELossBetheBlochTest.cxx | 12 ++++++------ AtTools/AtELossCATIMA.cxx | 11 +++++++---- scripts/formatAll.sh | 17 +++++++++++++---- 6 files changed, 43 insertions(+), 28 deletions(-) diff --git a/.github/workflows/CI-tests.yml b/.github/workflows/CI-tests.yml index 810be15de..33b8a9e01 100644 --- a/.github/workflows/CI-tests.yml +++ b/.github/workflows/CI-tests.yml @@ -18,7 +18,7 @@ jobs: runs-on: ubuntu-latest container: anthoak13/attpcroot-deps steps: - - uses: actions/checkout@v4 + - uses: actions/checkout@v6 - name: Configure CMake # Configure CMake in a 'build' subdirectory. `CMAKE_BUILD_TYPE` is only required if you are using a single-configuration generator such as make. diff --git a/.github/workflows/main.yml b/.github/workflows/main.yml index ffb510443..56253261b 100644 --- a/.github/workflows/main.yml +++ b/.github/workflows/main.yml @@ -8,7 +8,7 @@ jobs: runs-on: ubuntu-latest steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v6 - name: Run clang-format uses: DoozyX/clang-format-lint-action@v0.17 # These are optional (defaults displayed) diff --git a/AGENTS.md b/AGENTS.md index 1b5c6eb08..3b59135d5 100644 --- a/AGENTS.md +++ b/AGENTS.md @@ -4,21 +4,24 @@ ATTPCROOT is a ROOT/FairRoot-based C++ framework for simulation and analysis of ## Documentation -Full developer documentation lives in `docs/development/`: +Full developer documentation lives in `docs/`. See [docs/index.md](docs/index.md) for the full map. Quick topic links: | Topic | File | |-------|------| -| Environment setup | [docs/development/setup.md](../docs/development/setup.md) | -| Build & test commands | [docs/development/building.md](../docs/development/building.md) | -| Module overview | [docs/development/modules.md](../docs/development/modules.md) | -| Simulation pipeline | [docs/development/simulation-pipeline.md](../docs/development/simulation-pipeline.md) | -| Reconstruction pipeline | [docs/development/reconstruction-pipeline.md](../docs/development/reconstruction-pipeline.md) | -| Event generators | [docs/development/generators.md](../docs/development/generators.md) | -| Pulse shape analysis | [docs/development/psa.md](../docs/development/psa.md) | -| Energy loss | [docs/development/energy-loss.md](../docs/development/energy-loss.md) | -| Adding a new module | [docs/development/new-module.md](../docs/development/new-module.md) | -| Code style | [docs/development/code-style.md](../docs/development/code-style.md) | -| Testing | [docs/development/testing.md](../docs/development/testing.md) | +| First-time install | [tooling/installation.md](docs/tooling/installation.md) | +| Daily use (build/test) | [tooling/daily-use.md](docs/tooling/daily-use.md) | +| Testing patterns | [tooling/testing.md](docs/tooling/testing.md) | +| Contributor guide | [contributing/guide.md](docs/contributing/guide.md) | +| Adding a new module | [contributing/new-module.md](docs/contributing/new-module.md) | +| Code style | [contributing/code-style.md](docs/contributing/code-style.md) | +| Module overview | [reference/modules.md](docs/reference/modules.md) | +| Data model | [reference/data-model.md](docs/reference/data-model.md) | +| Branch I/O contracts | [reference/branch-io-contracts.md](docs/reference/branch-io-contracts.md) | +| Simulation pipeline | [subsystems/simulation-pipeline.md](docs/subsystems/simulation-pipeline.md) | +| Reconstruction pipeline | [subsystems/reconstruction-pipeline.md](docs/subsystems/reconstruction-pipeline.md) | +| Event generators | [subsystems/generators.md](docs/subsystems/generators.md) | +| Pulse shape analysis | [subsystems/psa.md](docs/subsystems/psa.md) | +| Energy loss | [subsystems/energy-loss.md](docs/subsystems/energy-loss.md) | ## Quick Reference: Build & Test diff --git a/AtTools/AtELossBetheBlochTest.cxx b/AtTools/AtELossBetheBlochTest.cxx index f1979fa75..f65ca93b3 100644 --- a/AtTools/AtELossBetheBlochTest.cxx +++ b/AtTools/AtELossBetheBlochTest.cxx @@ -1,4 +1,5 @@ #include "AtELossBetheBloch.h" + #include "AtELossCATIMA.h" #include @@ -286,17 +287,16 @@ TEST_F(AtELossBetheBlochVsCATIMAFixture, CATIMAComparison_EnergyStraggling) double eIni, eFin; }; const Case cases[] = { - {1.0, 0.75 * mass}, // narrow: ~25% energy loss - {5.0, 3.58164 * mass}, // narrow: ~28% energy loss - {5.0, 1.0}, // wide: 80% energy loss - {10.0, 1.0}, // wide: 90% energy loss, exercises Bragg-peak region + {1.0, 0.75 * mass}, // narrow: ~25% energy loss + {5.0, 3.58164 * mass}, // narrow: ~28% energy loss + {5.0, 1.0}, // wide: 80% energy loss + {10.0, 1.0}, // wide: 90% energy loss, exercises Bragg-peak region }; for (auto &c : cases) { double bbSigma = bb.GetElossStraggling(c.eIni, c.eFin); double catimaSigma = catima.GetElossStraggling(c.eIni, c.eFin); ASSERT_GT(bbSigma, 0.0) << "BB straggling zero for E0=" << c.eIni; ASSERT_GT(catimaSigma, 0.0) << "CATIMA straggling zero for E0=" << c.eIni; - EXPECT_NEAR(bbSigma, catimaSigma, 0.10 * catimaSigma) - << "Energy straggling mismatch for E0=" << c.eIni << " MeV"; + EXPECT_NEAR(bbSigma, catimaSigma, 0.10 * catimaSigma) << "Energy straggling mismatch for E0=" << c.eIni << " MeV"; } } diff --git a/AtTools/AtELossCATIMA.cxx b/AtTools/AtELossCATIMA.cxx index 2c26e5cd2..7920bdd58 100644 --- a/AtTools/AtELossCATIMA.cxx +++ b/AtTools/AtELossCATIMA.cxx @@ -44,7 +44,8 @@ double AtELossCATIMA::GetRange(double energyIni, double energyFin) const double remainingEnergy{energyIni}; double range{0}; while (remainingEnergy > energyFin) { - catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); + catima::Result result = + catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); double dEdx = result.dEdxi * fDensity; double DE = dEdx * fRangeStepSize / 10.; @@ -65,7 +66,8 @@ double AtELossCATIMA::GetEnergy(double energyIni, double distance) const double remainingEnergy{energyIni}; double range{0}; while (range < distance) { - catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); + catima::Result result = + catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); double dEdx = result.dEdxi * fDensity; double DE{}; @@ -142,8 +144,9 @@ AtELossCATIMA::GetBraggCurve(double energy, double rangeStepSize, double totalFr double range{}; while (remainingEnergy / energy > totalFractionELoss) { - catima::Result result = catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); - double dEdx = result.dEdxi * fDensity; // MeV/cm + catima::Result result = + catima::calculate(*fProjectile, *fMaterial, remainingEnergy / fProjectileMassAmu, fConfig); + double dEdx = result.dEdxi * fDensity; // MeV/cm braggCurve.push_back(std::make_pair(dEdx / 10.0, range)); // store in MeV/mm double DE = dEdx * rangeStepSize / 10.; diff --git a/scripts/formatAll.sh b/scripts/formatAll.sh index 873e68e9d..9a9594886 100755 --- a/scripts/formatAll.sh +++ b/scripts/formatAll.sh @@ -2,16 +2,25 @@ #Break if something errors set -e -echo "Indenting all source (.cxx) and header files (.h)" +echo "Indenting all source (.cxx, .C) and header files (.h)" -#Make sure $VMCWORKDIR is set. +#Make sure $VMCWORKDIR is set. if [ -z "$VMCWORKDIR" ] then echo "VMCWORKDIR is unset. Please run config.sh on build directory first. Aborting..." return fi -find -L $VMCWORKDIR -type f | grep -w h$ | xargs -I fname clang-format -i --verbose fname -find -L $VMCWORKDIR -type f | grep -w cxx$ | xargs -I fname clang-format -i --verbose fname +# Directories excluded from formatting (matches .github/workflows/main.yml) +EXCLUDE_DIRS="macro scripts compiled resources include geometry build .venv lib bin" + +# Build a prune expression for find +PRUNE_EXPR="" +for dir in $EXCLUDE_DIRS; do + PRUNE_EXPR="$PRUNE_EXPR -path $VMCWORKDIR/$dir -prune -o" +done + +find -L $VMCWORKDIR $PRUNE_EXPR -type f \( -name "*.h" -o -name "*.cxx" -o -name "*.C" \) -print \ + | xargs -I fname clang-format-17 -i --verbose fname echo "Success! Indented all source and header files" From eb49be7ae3c608967fab67f6ca20b86590afa1de Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Wed, 1 Apr 2026 18:23:38 -0400 Subject: [PATCH 07/10] Update documentation --- docs/contributing/code-style.md | 33 +++++- docs/index.md | 6 + docs/reference/branch-io-contracts.md | 49 ++++---- docs/reference/data-model.md | 46 ++++--- docs/reference/macro-cookbook.md | 2 +- docs/reference/modules.md | 22 ++-- docs/reference/parameters.md | 93 +++++++++++++++ docs/subsystems/geometry.md | 132 +++++++++++++++++++++ docs/subsystems/psa.md | 24 ++-- docs/subsystems/reconstruction-pipeline.md | 21 ++-- docs/subsystems/simulation-pipeline.md | 77 ++++++++---- docs/subsystems/visualization.md | 112 +++++++++++++++++ docs/tooling/testing.md | 21 ++++ 13 files changed, 538 insertions(+), 100 deletions(-) create mode 100644 docs/reference/parameters.md create mode 100644 docs/subsystems/geometry.md create mode 100644 docs/subsystems/visualization.md diff --git a/docs/contributing/code-style.md b/docs/contributing/code-style.md index 33438df3a..e0062ee0b 100644 --- a/docs/contributing/code-style.md +++ b/docs/contributing/code-style.md @@ -1,6 +1,37 @@ # Code Style -This page covers formatting and static-analysis tools only. Contributor-critical ROOT and persistence rules live in [guide.md](guide.md). +This page covers naming conventions, formatting, and static-analysis tools. Contributor-critical ROOT and persistence rules live in [guide.md](guide.md). + +## Naming Conventions + +### Libraries and Classes + +- **Libraries**: a folder named `AtFoo` produces shared object `libAtFoo` +- **Classes**: non-namespaced classes use the `At` prefix in CamelCase — `AtBaseEvent`, `AtPSAMax`, `AtFitterTask` +- Namespaced code (e.g., `EventFit::`, `kf::`) omits the `At` prefix + +### Data Members + +- Private and protected members: `f` prefix + capital letter — `fTrackID`, `fDriftVelocity` +- Boolean flags: `k` prefix — `kIsGood`, `kInitialized` +- Local variables and function parameters: no prefix, lowercase camelCase + +### Methods + +Follow ROOT naming conventions: +- Getters: `GetXaxis()` not `GetXAxis()` (avoid back-to-back capitals) +- Setters: `SetThreshold()` +- Boolean queries: `IsValid()`, `HasHits()` + +### ClassDef / ClassImp + +Use `ClassDef` and `ClassImp` **only** for `TObject` subclasses. When the memory layout changes, increment the version number in `ClassDef`: + +```cpp +ClassDef(AtMyData, 2); // bump when adding/removing/reordering members +``` + +Non-persisted classes (tasks, algorithms, models) do not need `ClassDef`. ## Formatting diff --git a/docs/index.md b/docs/index.md index 2a6de11ac..e8f9c43ab 100644 --- a/docs/index.md +++ b/docs/index.md @@ -13,6 +13,9 @@ Agent-facing docs for the current branch state of ATTPCROOT. - understand runtime objects: [reference/data-model.md](reference/data-model.md) - find subsystem ownership: [reference/modules.md](reference/modules.md) - find an example macro: [reference/macro-cookbook.md](reference/macro-cookbook.md) +- configure detector parameters: [reference/parameters.md](reference/parameters.md) +- prepare detector geometry: [subsystems/geometry.md](subsystems/geometry.md) +- browse data interactively: [subsystems/visualization.md](subsystems/visualization.md) ## Architecture @@ -23,6 +26,9 @@ Agent-facing docs for the current branch state of ATTPCROOT. - generators: [subsystems/generators.md](subsystems/generators.md) - PSA: [subsystems/psa.md](subsystems/psa.md) - energy loss: [subsystems/energy-loss.md](subsystems/energy-loss.md) +- geometry: [subsystems/geometry.md](subsystems/geometry.md) +- visualization: [subsystems/visualization.md](subsystems/visualization.md) +- parameters: [reference/parameters.md](reference/parameters.md) ## Branch-Specific Notes diff --git a/docs/reference/branch-io-contracts.md b/docs/reference/branch-io-contracts.md index a1d69c09c..71d1e5d97 100644 --- a/docs/reference/branch-io-contracts.md +++ b/docs/reference/branch-io-contracts.md @@ -1,33 +1,40 @@ # Branch and IO Contracts -This page records the current branch names and persistence defaults used by major tasks. These defaults come from the task constructors and `Init()` implementations in the current branch. +Default branch names, contained types, and persistence settings for major tasks. These come from task constructors and `Init()` implementations in the current source. -In most cases the branch object returned by `FairRootManager::GetObject(...)` is a `TClonesArray`, not a bare event object. The task then reads slot `0` as the current `AtRawEvent`, `AtEvent`, `AtPatternEvent`, or `AtTrackingEvent`. +Every branch listed here is a `TClonesArray` holding one event object at slot `0` — not a bare pointer: + +```cpp +auto *arr = dynamic_cast(ioMan->GetObject("AtEventH")); +auto *event = dynamic_cast(arr->At(0)); +``` ## Reconstruction Tasks -| Task | Input Branches | Output Branches | Persistence Default | Notes | -|------|----------------|-----------------|---------------------|-------| -| `AtUnpackTask` | external file via unpacker | `AtRawEvent` | `true` | Registers a `TClonesArray` holding `AtRawEvent` objects | -| `AtFilterTask` | `AtRawEvent` | `AtRawEventFiltered` | `false` | Reads/writes `TClonesArray`-style branches; can also filter aux/FPN pads | -| `AtPSAtask` | `AtRawEvent` | `AtEventH` | `false` | Reads `TClonesArray` input, writes `TClonesArray("AtEvent", 1)` output; optionally looks for `AtTpcPoint` | -| `AtDataCleaningTask` | `AtEventH` | `AtEventCleaned` | `true` | Reads the input event from a `TClonesArray` and writes a cleaned `AtEvent` copy | -| `AtPRAtask` | `AtEventH` | `AtPatternEvent` | `false` | PRA-based pattern recognition path using `TClonesArray` event containers | -| `AtSampleConsensusTask` | `AtEventH` | `AtPatternEvent` | `false` | Sample-consensus path using the same container pattern | -| `AtFitterTask` | `AtPatternEvent` | `AtTrackingEvent` | `false` | Output is a `TClonesArray("AtTrackingEvent", 1)`; verify input-branch behavior in source before relying on setters | -| `AtMCFitterTask` | `AtPatternEvent` | `AtMCResult`, `SimEvent`, `SimRawEvent` | `true`, `false`, `false` | Save toggles exist for each output branch | +| Task | Input branch / type | Output branch / type | Persist default | Branch setters | +|------|---------------------|----------------------|-----------------|----------------| +| `AtUnpackTask` | external file via unpacker | `AtRawEvent` → `TClonesArray` | **true** | — | +| `AtFilterTask` | `AtRawEvent` → `TClonesArray` | `AtRawEventFiltered` → `TClonesArray` | false | `SetInputBranch`, `SetOutputBranch` | +| `AtPSAtask` | `AtRawEvent` → `TClonesArray` | `AtEventH` → `TClonesArray` | false | `SetInputBranch`, `SetOutputBranch` | +| `AtDataCleaningTask` | `AtEventH` → `TClonesArray` | `AtEventCleaned` → `TClonesArray` | **true** | `SetInputBranch`, `SetOutputBranch` | +| `AtPRAtask` | `AtEventH` → `TClonesArray` | `AtPatternEvent` → `TClonesArray` | false | `SetInputBranch`, `SetOutputBranch` | +| `AtSampleConsensusTask` | `AtEventH` → `TClonesArray` | `AtPatternEvent` → `TClonesArray` | false | `SetInputBranch`, `SetOutputBranch` | +| `AtFitterTask` | `AtPatternEvent` → `TClonesArray` | `AtTrackingEvent` → `TClonesArray` | false | `SetInputBranch`, `SetOutputBranch` | +| `AtMCFitterTask` | `AtPatternEvent` → `TClonesArray` | `AtMCResult` (**true**), `SimEvent` (false), `SimRawEvent` (false) | per-branch | `SetPatternBranchName`, `SetSaveResult/Event/RawEvent` | ## Simulation and Digitization Tasks -| Task | Input Branches | Output Branches | Persistence Default | Notes | -|------|----------------|-----------------|---------------------|-------| -| `AtClusterizeTask` | `AtTpcPoint` | `AtSimulatedPoint` | `false` | Reads the detector-specific `AtTpcPoint` `TClonesArray`, not a generic `AtMCPoint` branch | -| `AtPulseTask` | `AtSimulatedPoint`, optional `AtTpcPoint` | `AtRawEvent`, optional re-registered `AtTpcPoint` | `true`, `false` | Reads branch containers, not bare objects; `SetSaveMCInfo()` enables MC mapping behavior | +| Task | Input branch / type | Output branch / type | Persist default | Branch setters | +|------|---------------------|----------------------|-----------------|----------------| +| `AtClusterizeTask` | `AtTpcPoint` → `TClonesArray` | `AtSimulatedPoint` → `TClonesArray` | false | `SetPersistence` | +| `AtPulseTask` | `AtSimulatedPoint` → `TClonesArray`, optional `AtTpcPoint` → `TClonesArray` | `AtRawEvent` → `TClonesArray` (**true**); re-registers `AtTpcPoint` (false) | true / false | `SetOutputBranch`, `SetPersistence`, `SetPersistenceAtTpcPoint`, `SetSaveMCInfo` | ## Reading This Page -- These are default names, not hard-coded universal laws. -- Many tasks expose `SetInputBranch(...)`, `SetOutputBranch(...)`, or related setters. -- Do not assume every setter is honored consistently; check the task source before automating branch-name rewrites. -- If you add a task that reads/writes FairRoot branches, update this page. -- For the runtime objects behind these branches, see [data-model.md](data-model.md). +- These are defaults, not hard-coded universals. Most tasks expose `SetInputBranch` / `SetOutputBranch`. +- `AtPRAtask` and `AtSampleConsensusTask` write the same output branch name — only one should be active per run. +- `AtMCFitterTask` writes three separate branches; each has its own save toggle. +- `AtFittedTrack` (inside `AtTrackingEvent`) is currently unstable — treat layout as best-effort. +- In the AT-TPC path, `AtTpcPoint` is a branch name whose current container holds `AtMCPoint` objects. +- If you add a task that reads or writes FairRoot branches, add a row above. +- For what the contained objects hold (fields, ownership), see [data-model.md](data-model.md). diff --git a/docs/reference/data-model.md b/docs/reference/data-model.md index a40a1fb6a..2b5366231 100644 --- a/docs/reference/data-model.md +++ b/docs/reference/data-model.md @@ -6,41 +6,37 @@ Most FairRoot branches in this tree are not bare `AtEvent*` or `AtRawEvent*` obj ## Reconstruction Objects -| Object | Typical branch container | Produced by | Consumed by | Persisted | Key downstream fields | -|--------|--------------------------|-------------|-------------|-----------|-----------------------| -| `AtRawEvent` | `TClonesArray` branch such as `AtRawEvent` / `AtRawEventFiltered` | `AtUnpackTask`, `AtPulseTask` | filters, PSA | yes | pad traces, aux/FPN pads, good flag, optional MC map | -| `AtEvent` | `TClonesArray` branch such as `AtEventH` / `AtEventCleaned` | `AtPSAtask` | cleaning, pattern recognition | yes | hit list, event charge, mesh signal | -| `AtPatternEvent` | `TClonesArray` branch `AtPatternEvent` | `AtPRAtask`, `AtSampleConsensusTask` | fitting | yes | candidate tracks, noise hits | -| `AtTrack` | contained inside `AtPatternEvent`, not its own top-level branch | pattern recognition code | fitters | yes | track ID, hit/cluster-hit collections, geometric estimates, Bragg-curve values | +| Object | Typical branch container | Main producer | Main consumer | Persisted | Key downstream fields | +|--------|--------------------------|---------------|---------------|-----------|-----------------------| +| `AtRawEvent` | `TClonesArray` branch such as `AtRawEvent` / `AtRawEventFiltered` | `AtUnpackTask`, `AtPulseTask` | `AtFilterTask`, `AtPSAtask` | yes | pad traces, aux/FPN pads, good flag, optional MC map | +| `AtEvent` | `TClonesArray` branch such as `AtEventH` / `AtEventCleaned` | `AtPSAtask` | `AtDataCleaningTask`, `AtPRAtask`, `AtSampleConsensusTask` | yes | hit list, event charge, mesh signal | +| `AtPatternEvent` | `TClonesArray` branch such as `AtPatternEvent` / `AtPatternEventModified` | `AtPRAtask`, `AtSampleConsensusTask` | `AtFitterTask`, `AtMCFitterTask` | yes | candidate tracks, noise hits | +| `AtTrack` | contained inside `AtPatternEvent`, not its own top-level branch | pattern-recognition code | fitters, track-level analysis | yes | track ID, hit/cluster-hit collections, geometric estimates, Bragg-curve values | | `AtTrackingEvent` | `TClonesArray` branch `AtTrackingEvent` | `AtFitterTask` | downstream analysis | yes | fitted tracks, optional copied track array, event vertex fields | | `AtFittedTrack` | contained inside `AtTrackingEvent`, not its own top-level branch | fitter implementations | downstream analysis | yes | current branch layout only; treat as unstable | +| `AtMCResult` | `TClonesArray` branch `AtMCResult` | `AtMCFitterTask` | MC-fitting analysis | yes | objective value and fitted/sampled parameter summary | ## Simulation Objects -| Object | Typical branch container | Produced by | Consumed by | Persisted | Notes | -|--------|--------------------------|-------------|-------------|-----------|-------| -| `AtTpcPoint` | `TClonesArray` branch `AtTpcPoint` | Geant/VMC transport through the AT-TPC detector | digitization tasks such as `AtClusterizeTask` and `AtPulseTask` | yes | detector-specific point branch seen by downstream digitization code | -| `AtMCPoint` | module-specific MC point branches | Geant/VMC transport | analysis, truth matching, detector-specific point types | yes | generic MC-point base family; downstream digitization usually works with detector-specific derived point classes | +| Object | Typical branch container | Main producer | Main consumer | Persisted | Notes | +|--------|--------------------------|---------------|---------------|-----------|-------| +| `AtMCPoint` | MC-point branches such as `AtTpcPoint` | Geant/VMC transport, `AtSimpleSimulation` | `AtClusterizeTask` | yes | primary simulation-to-digitization handoff object; also used for MC bookkeeping in simulated-data reconstruction | +| `AtSimulatedPoint` | `TClonesArray` branch `AtSimulatedPoint` | `AtClusterizeTask` | `AtPulseTask` | yes | intermediate ionization-electron / charge-cluster representation between MC points and pad traces | | `AtMCTrack` | `TClonesArray` branch `MCTrack` | simulation stack | analysis, truth matching | yes | simulated particle tracks | | `AtVertexPropagator` | no FairRoot branch | generators/runtime simulation logic | generators and downstream simulation code | no | singleton shared state | -## Flow Summary +## Container Pattern ```text -TClonesArray(AtRawEvent branch) -> AtRawEvent -TClonesArray(AtEventH branch) -> AtEvent -TClonesArray(AtPatternEvent) -> AtPatternEvent -TClonesArray(AtTrackingEvent) -> AtTrackingEvent - | | - | -> contains AtTrack - -> contains hits -``` - -```text -Generator / AtVertexPropagator - -> MCTrack / AtTpcPoint branches - -> AtSimulatedPoint branch - -> AtRawEvent branch +FairRoot branch name -> TClonesArray container -> runtime object + +AtRawEvent -> TClonesArray -> AtRawEvent +AtEventH -> TClonesArray -> AtEvent +AtPatternEvent -> TClonesArray -> AtPatternEvent +AtTrackingEvent -> TClonesArray -> AtTrackingEvent +MCTrack -> TClonesArray -> AtMCTrack +AtTpcPoint -> TClonesArray -> AtMCPoint +AtSimulatedPoint-> TClonesArray -> AtSimulatedPoint ``` For the tasks that produce and consume these objects, see [branch-io-contracts.md](branch-io-contracts.md). For pipeline order, see [simulation-pipeline.md](../subsystems/simulation-pipeline.md) and [reconstruction-pipeline.md](../subsystems/reconstruction-pipeline.md). diff --git a/docs/reference/macro-cookbook.md b/docs/reference/macro-cookbook.md index 0ff481383..6e4938fba 100644 --- a/docs/reference/macro-cookbook.md +++ b/docs/reference/macro-cookbook.md @@ -1,6 +1,6 @@ # Macro Cookbook -Example macro starting points referenced by these docs. This page is intentionally sparse. +Example macro starting points referenced by these docs. ## Notes diff --git a/docs/reference/modules.md b/docs/reference/modules.md index 4182e3228..919e81a0f 100644 --- a/docs/reference/modules.md +++ b/docs/reference/modules.md @@ -4,18 +4,20 @@ ATTPCROOT is organized as CMake library targets layered around data, simulation, ## Modules +`Also inspect` points to the most useful overview docs for that module. + | Module | Role | Main subareas | Common edit reasons | Also inspect | |--------|------|---------------|---------------------|--------------| -| `AtData` | persisted event and track containers | event, track, pattern classes | add/modify data objects | LinkDef, producing tasks, [data-model.md](data-model.md) | -| `AtReconstruction` | PSA, filtering, pattern recognition, fitting tasks | `AtPulseAnalyzer/`, `AtPatternRecognition/`, `AtFitter/`, `AtFilter/` | change branch flow or task behavior | [branch-io-contracts.md](branch-io-contracts.md), tests | -| `AtUnpack` | experimental data to `AtRawEvent` | unpackers, `GETDecoder2/` | add file-format/front-end support | maps, macros | -| `AtDigitization` | simulation truth to detector signals | clusterization, pulse generation | change MC-to-signal behavior | simulation pipeline | -| `AtTools` | shared utilities | energy loss, kinematics, hit sampling, cleaning | algorithmic utilities used across modules | dependent modules | -| `AtSimulationData` | MC truth and simulation shared state | MC objects, propagator | change truth object layout or simulation state | simulation pipeline | -| `AtGenerators` | FairRoot generator implementations | reaction and beam generators | change simulation source behavior | `AtVertexPropagator`, simulation pipeline | -| `AtDetectors` | detector geometry/sensitive detectors | detector-specific subdirs | detector implementation changes | geometry assets | -| `AtAnalysis` | higher-level analysis code | experiment-facing analysis | analysis task changes | data model | -| `AtEventDisplay` | event visualization | tabs/display code | visualization changes | data model | +| `AtData` | persisted event and track containers | event, track, pattern classes | add/modify data objects | [data-model.md](data-model.md), [branch-io-contracts.md](branch-io-contracts.md) | +| `AtReconstruction` | PSA, filtering, pattern recognition, fitting tasks | `AtPulseAnalyzer/`, `AtPatternRecognition/`, `AtPatternModification/`, `AtFitter/`, `AtFilter/` | change branch flow or task behavior | [reconstruction-pipeline.md](../subsystems/reconstruction-pipeline.md), [branch-io-contracts.md](branch-io-contracts.md) | +| `AtUnpack` | experimental data to `AtRawEvent` | unpackers, `GETDecoder2/`, `deprecated/` | add file-format/front-end support | [data-model.md](data-model.md) | +| `AtDigitization` | simulation truth to detector signals | clusterization, pulse generation, trigger/space-charge tasks | change MC-to-signal behavior | [simulation-pipeline.md](../subsystems/simulation-pipeline.md), [branch-io-contracts.md](branch-io-contracts.md) | +| `AtTools` | shared utilities | energy loss, kinematics, hit sampling, cleaning, electronic response | algorithmic utilities used across modules | [energy-loss.md](../subsystems/energy-loss.md) | +| `AtSimulationData` | MC truth and simulation shared state | MC objects, stack, propagator | change truth object layout or simulation state | [simulation-pipeline.md](../subsystems/simulation-pipeline.md), [data-model.md](data-model.md) | +| `AtGenerators` | FairRoot generator implementations | reaction and beam generators | change simulation source behavior | [generators.md](../subsystems/generators.md), [simulation-pipeline.md](../subsystems/simulation-pipeline.md) | +| `AtDetectors` | detector geometry and sensitive detectors | detector-specific subdirs, field/passive geometry | detector implementation changes | [geometry.md](../subsystems/geometry.md) | +| `AtAnalysis` | higher-level analysis code | experiment-facing analysis | analysis task changes | [data-model.md](data-model.md), [reconstruction-pipeline.md](../subsystems/reconstruction-pipeline.md) | +| `AtEventDisplay` | event visualization | `AtTabs/`, `AtSidebar/`, legacy display code | visualization changes | [visualization.md](../subsystems/visualization.md), [data-model.md](data-model.md) | ## Workflow Assets diff --git a/docs/reference/parameters.md b/docs/reference/parameters.md new file mode 100644 index 000000000..cacd18d88 --- /dev/null +++ b/docs/reference/parameters.md @@ -0,0 +1,93 @@ +# Parameter Files + +ATTPCROOT uses FairRoot's parameter system to manage experiment-specific detector configuration. The main container class is `AtDigiPar`, defined in `AtParameter/AtDigiPar.h`. + +## Loading Parameters in a Macro + +```cpp +FairRuntimeDb *rtdb = run->GetRuntimeDb(); +FairParAsciiFileIo *parIo = new FairParAsciiFileIo(); +parIo->open("parameters/ATTPC.myexperiment.par", "in"); +rtdb->setFirstInput(parIo); +rtdb->getContainer("AtDigiPar"); +``` + +Parameter files live in `parameters/`. Each file typically contains an `[AtDigiPar]` section followed by key-value pairs. Pick the file that matches the experiment or create a new one by copying an existing one. + +## AtDigiPar Parameters + +### Electromagnetic Fields + +| Parameter | Type | Unit | Description | +|-----------|------|------|-------------| +| `EField` | Double_t | V/m | Longitudinal drift electric field | +| `BField` | Double_t | T | Longitudinal magnetic field | + +### Detector Geometry + +| Parameter | Type | Unit | Description | +|-----------|------|------|-------------| +| `TBEntrance` | Int_t | time buckets | Beam entrance position at detector entrance | +| `ZPadPlane` | Double_t | mm | Position of the micromegas pad plane (detector length) | +| `DriftLength` | Double_t | mm | Full drift length (AT-TPC: 1000, prototype: 500) | + +### Gas Physics + +| Parameter | Type | Unit | Description | +|-----------|------|------|-------------| +| `DriftVelocity` | Double_t | cm/µs | Electron drift velocity | +| `EIonize` | Double_t | eV | Effective ionization energy of the fill gas | +| `Fano` | Double_t | — | Fano factor of the fill gas | +| `CoefL` | Double_t | cm^−½ | Longitudinal diffusion coefficient | +| `CoefT` | Double_t | cm^−½ | Transverse diffusion coefficient | +| `GasPressure` | Double_t | torr | Gas pressure | +| `Density` | Double_t | kg/m³ | Gas density | + +### Electronics + +| Parameter | Type | Unit | Description | +|-----------|------|------|-------------| +| `SamplingRate` | Int_t | MHz | GET electronics sampling frequency | +| `GETGain` | Double_t | fC | Gain from GET electronics | +| `PeakingTime` | Int_t | ns | Electronic response peaking time | +| `Gain` | Double_t | — | Average micromegas amplification factor | +| `NumTbs` | Int_t | — | Number of time buckets | + +## Accessing Parameters in Code + +Tasks retrieve `AtDigiPar` from `FairRuntimeDb` in their `Init()`: + +```cpp +auto *rtdb = FairRuntimeDb::instance(); +fPar = dynamic_cast(rtdb->getContainer("AtDigiPar")); +``` + +Then in `Exec()`: + +```cpp +double driftVel = fPar->GetDriftVelocity(); // cm/µs +double bField = fPar->GetBField(); // T +``` + +## Example Parameter File Entry + +``` +[AtDigiPar] +EField:Double_t 5000 # V/m +BField:Double_t 2 # Tesla +TBEntrance:Int_t 280 # time bucket at detector entrance +ZPadPlane:Double_t 1000 # mm +DriftVelocity:Double_t 5.00 # cm/us +EIonize:Double_t 15.603 # eV (He+CO2 example) +Fano:Double_t 0.22 +CoefL:Double_t 0.025 # cm^-0.5 +CoefT:Double_t 0.010 # cm^-0.5 +GasPressure:Double_t 100 # torr +Density:Double_t 0.0738 # kg/m3 +SamplingRate:Int_t 12.5 # MHz +Gain:Double_t 100.0 +NumTbs:Int_t 512 +PeakingTime:Int_t 117 # ns +``` + +See `parameters/AT.parameters.par` for a complete reference example and `parameters/` for all experiment-specific files. diff --git a/docs/subsystems/geometry.md b/docs/subsystems/geometry.md new file mode 100644 index 000000000..a52ee6581 --- /dev/null +++ b/docs/subsystems/geometry.md @@ -0,0 +1,132 @@ +# Detector Geometry + +Use this file when the task involves detector geometry assets, simulation geometry loading, or geometry-related navigation. + +## Start Here + +- `AtDetectors/AtTpc/AtTpc.cxx` +- `AtDetectors/AtTpc/AtTpc.h` +- `geometry/` +- `macro/examples/run_sim.C` +- `macro/tests/AT-TPC/run_sim_attpc.C` +- `macro/tests/generateGeometry.sh` + +If the task is about pad centers, pad IDs, or pad-plane layout, also inspect: + +- `AtMap/AtTpcMap.cxx` +- `AtMap/AtSpecMATMap.cxx` + +## Runtime Path + +Simulation macros pass a geometry file to `AtTpc` with `SetGeometryFileName()`: + +```cpp +AtTpc *tpc = new AtTpc("ATTPC", kTRUE); +tpc->SetGeometryFileName("ATTPC_v1.2.root"); +run->AddModule(tpc); +``` + +`AtTpc::ConstructGeometry()` dispatches by extension: + +- `.root` -> `ConstructRootGeometry()` +- `.geo` -> ASCII branch exists but is not the practical path here + +For this tree, `AtTpc` geometry is normally loaded from generated ROOT files. + +## Where Geometry Lives + +`geometry/` contains: + +- generator macros such as `ATTPC_v1_2.C`, `ATTPC_He1bar.C`, `GADGET_II.C` +- generated outputs such as `ATTPC_v1.2.root` and `ATTPC_v1.2_geomanager.root` +- `media.geo` for material definitions +- detector-specific helper assets + +Do not assume geometry outputs are missing. This checkout already contains several generated `.root` geometry files. + +## Generation + +Generate a geometry only if the required output is missing or stale. + +Example: + +```bash +cd geometry +root -l -q ATTPC_v1_2.C +``` + +This writes: + +- `ATTPC_v1.2.root` +- `ATTPC_v1.2_geomanager.root` + +Test helper: + +```bash +macro/tests/generateGeometry.sh +``` + +That script generates only a fixed subset of geometries used by tests. It is not a full geometry build step. + +## Detector Variants + +| Script | Detector | Notes | +|--------|----------|-------| +| `ATTPC_v1_2.C` | Full AT-TPC | Used by most unpacking/display macros | +| `ATTPC_v1_1.C` | Full AT-TPC | Older version; some example sim macros still reference this | +| `ATTPC_Proto_v1_0.C` | Prototype AT-TPC | Shorter drift length (500 mm) | +| `ATTPC_He1bar.C` | AT-TPC, He 1 atm | Used by the AT-TPC integration test | +| `ATTPC_He*.C`, `ATTPC_D*.C`, `ATTPC_H*.C` | Gas/pressure variants | Specific gas fill and pressure combinations | +| `SpecMAT_*.C` | SpecMAT | Different pad geometry; scintillator variants available | +| `GADGET_II.C` | GADGET-II | Cylindrical TPC variant | +| `FRIB_TPC_v*.C` | FRIB TPC | | + +## Macro Anchors + +Use real macros, not isolated snippets, when tracing behavior: + +- `macro/examples/run_sim.C` +- `macro/tests/AT-TPC/run_sim_attpc.C` + +Those show the surrounding `FairRunSim` setup, `run->SetMaterials("media.geo")`, and the other modules commonly registered with the TPC geometry such as `AtCave`, `AtMagnet`, and `AtPipe`. + +## Variant Selection + +There is no single geometry file that is clearly canonical across the tree. + +Different workflows use different files: + +- unpacking and display macros often use `ATTPC_v1.2.root` +- older example simulation macros use `ATTPC_v1.1.root` +- the AT-TPC integration test uses `ATTPC_He1bar.root` + +When editing a workflow, trust the specific macro or test in that workflow over generic prose. + +## Materials + +`run->SetMaterials("media.geo")` supplies material definitions used by geometry generation and simulation setup. + +Common modifications to `geometry/media.geo`: + +- **Gas composition**: update the `AtTPC_Gas` entry — change density, atomic composition, and weight fractions +- **Window materials**: modify the entrance window material entry +- **Add a material**: append a new block following the existing format; the file header documents the required fields (`ncomp`, `aw`, `an`, `dens`, `wm`, `sensflag`, `fldflag`, `fld`, `epsil`) + +If `media.geo` changes, regenerate any geometry outputs that depend on it. + +## Geometry vs Mapping + +Keep these separate: + +- transport geometry: `AtDetectors/AtTpc`, `geometry/*.C`, simulation macros +- pad mapping and display-plane geometry: `AtMap/`, visualization, decoder/mapping setup + +If the task mentions pad geometry or lookup tables, start in `AtMap/`, not `AtDetectors/AtTpc`. + +## Navigation Order + +1. Find the macro, task, or test that names the geometry file. +2. Check whether that `.root` file already exists in `geometry/`. +3. Inspect `AtDetectors/AtTpc/AtTpc.cxx` for the runtime load path. +4. Inspect the matching geometry generator macro in `geometry/` if the asset must change. +5. Branch into `AtMap/` only if the task is really about pad mapping rather than transport geometry. diff --git a/docs/subsystems/psa.md b/docs/subsystems/psa.md index 2cfae634a..13af6e209 100644 --- a/docs/subsystems/psa.md +++ b/docs/subsystems/psa.md @@ -6,11 +6,17 @@ PSA extracts hit positions and charge information from raw pad traces (`AtRawEve `AtPSA` is the abstract base class in `AtReconstruction/AtPulseAnalyzer/`. All concrete PSA implementations inherit from it. -| Class | Algorithm | -|-------|-----------| -| `AtPSAMax` | Peak-finding: hit at the time bucket of maximum signal | -| `AtPSAFull` | Integrates the full trace; produces one hit per pad | -| `AtPSADeconv` | Deconvolution-based hit finding; can find multiple hits per pad | +| Class | Hits per pad | Description | +|-------|-------------|-------------| +| `AtPSAMax` | 1 | Places one hit at the time bucket of peak ADC. Charge stored is the peak ADC value; total trace integral also stored on the hit. Optional sub-TB time correction via charge-weighted average of peak region. | +| `AtPSAFull` | 1 or many | Designed for long traces (e.g. beam tracks). Traces with a signal window <50 TBs get one hit at the peak; traces ≥50 TBs are split into 50-TB windows with one hit per window and charge = average ADC in that window. | +| `AtPSASpectrum` | 1 or many | Uses ROOT `TSpectrum` for peak finding; can find multiple peaks per pad. Supports optional background subtraction and interpolation. | +| `AtPSAHitPerTB` | many | One hit per time bucket above threshold. Produces a dense 3D charge distribution with no peak-finding assumptions. Configurable TB range. | +| `AtPSATBAvg` | many | Divides the trace into windows of N TBs (default 5) and creates one hit per window whose average charge exceeds threshold. Hit Z is at window center. Skips saturated pads when a max threshold is set. Can operate on a pad augment instead of the raw ADC. | +| `AtPSADeconv` | 1 | FFT-deconvolves the electronics response function from the raw trace and applies a Butterworth low-pass filter to recover the ionization current. Reconstructed charge saved as `"Qreco"` augment on the pad. Hit placed at charge-distribution weighted average. Requires a response function (`SetResponse`). | +| `AtPSADeconvFit` | 1 | Extends `AtPSADeconv`; instead of the weighted average, fits the deconvolved charge distribution with a Gaussian using the longitudinal diffusion coefficient for improved Z resolution. | +| `AtPSAIterDeconv` | 1 | Extends `AtPSADeconv` with iterative corrections to the reconstructed current. | +| `AtPSASimple2` | — | **Deprecated.** Use `AtPSASpectrum` or `AtPSAMax` instead. | ## Running PSA @@ -18,7 +24,10 @@ PSA is applied by `AtPSAtask`, which wraps any `AtPSA` subclass: ```cpp auto psa = std::make_unique(); +psa->SetThreshold(45); // ADC threshold; optional auto psaTask = new AtPSAtask(std::move(psa)); +psaTask->SetInputBranch("AtRawEventFiltered"); // default: "AtRawEvent" +psaTask->SetPersistence(true); fRun->AddTask(psaTask); ``` @@ -36,5 +45,6 @@ composite->AddPSA(std::make_unique()); PSA produces an `AtEvent` stored as a branch in the output TTree. Each `AtHit` in the event carries: - 3D position (x, y, time-bucket → z via drift velocity) -- Integrated charge -- Pad number and other metadata +- Charge (meaning depends on algorithm — peak ADC, window average, etc.) +- Trace integral (sum of all ADC values over the trace) +- Pad number and timestamp metadata diff --git a/docs/subsystems/reconstruction-pipeline.md b/docs/subsystems/reconstruction-pipeline.md index 4845df473..450286968 100644 --- a/docs/subsystems/reconstruction-pipeline.md +++ b/docs/subsystems/reconstruction-pipeline.md @@ -28,21 +28,16 @@ The exact branch names and defaults are documented in [branch-io-contracts.md](. ## Pattern Recognition -Two supported paths currently produce `AtPatternEvent`: +Two supported paths produce `AtPatternEvent → TClonesArray[AtPatternEvent]`. Both populate `fHitArray` on each track and `fNoise` on the event with unassigned hits. They differ in which `AtTrack` fields they fill: -```text -AtPRAtask - │ PRA / triplclust-style pattern recognition - └─ output branch: AtPatternEvent -> TClonesArray[AtPatternEvent] -``` - -```text -AtSampleConsensusTask - │ sample-consensus pattern recognition - └─ output branch: AtPatternEvent -> TClonesArray[AtPatternEvent] -``` +| | `AtPRAtask` | `AtSampleConsensusTask` | +|-|-------------|-------------------------| +| Algorithm | triplet/hierarchical clustering (TriplClust default; `SetPRAlgorithm(n)`) | sample consensus (`SetPatternType`, `SetEstimator`) | +| `fPattern` | always `AtPatternCircle2D` | configurable: `kLine`, `kCircle2D`, `kRay`, `kY` | +| `fGeoRadius`, `fGeoCenter`, `fGeoTheta`, `fGeoPhi` | **filled** — via internal RANSAC after clustering | **not filled** (remain zero) | +| `fHitClusterArray` | populated | not populated | -These are both supported. The docs do not declare one of them to be the universal default. +The `fGeo*` fields are the momentum-seed inputs for GenFit-based fitters. `AtSampleConsensusTask` leaves them empty, so it requires additional processing (e.g. `AtPatternModificationTask`) before a GenFit-based fitter can run. ## Fitting diff --git a/docs/subsystems/simulation-pipeline.md b/docs/subsystems/simulation-pipeline.md index cccb33414..e808e922c 100644 --- a/docs/subsystems/simulation-pipeline.md +++ b/docs/subsystems/simulation-pipeline.md @@ -2,53 +2,86 @@ The simulation pipeline converts a generated reaction into MC truth, simulated detector responses, and finally `AtRawEvent` output that can be passed into reconstruction. -Like reconstruction, the task-facing branch objects here are usually `TClonesArray` containers. The detector-point branch used by digitization is `AtTpcPoint`, not a generic `AtMCPoint` branch. +There are currently two ways to generate the simulation-side `AtMCPoint` input: -## Steps +- the full FairRoot/VMC transport path used by most detector simulations +- a lighter-weight model-driven path built around `AtSimpleSimulation` + +## Flow + +The overall simulation flow is shared downstream of `AtMCPoint` generation: ``` -FairPrimaryGenerator - └─ AtReactionGenerator subclass - │ generates primary particles and writes to particle stack - ▼ -Geant4 / VMC transport - │ produces detector-point branches such as AtTpcPoint +FairPrimaryGenerator + AtReactionGenerator AtSimpleSimulation + │ │ + ▼ ▼ + Geant4 / VMC transport direct AtMCPoint generation + │ │ + └───────────────┬──────────────────┘ + ▼ + AtMCPoint ▼ AtClusterizeTask - │ reads TClonesArray(AtTpcPoint) and converts deposits -> ionization electron clusters + │ converts MC-point deposits -> ionization electron clusters ▼ AtPulseTask │ drifts electrons and produces simulated pad traces └─ output branch: AtRawEvent -> TClonesArray[AtRawEvent] ``` +The only difference between the two paths is how `AtMCPoint` objects are generated. After that, the downstream digitization flow is the same. + +### AtMCPoint Generation + +- **Primary path:** Geant4/VMC transport through the detector geometry produces `AtMCPoint` objects from the generated particles. +- **Secondary path:** `AtDigitization/AtSimpleSimulation.h` advances particles through the active volume in fixed steps, applies an `AtTools::AtELossModel`, and writes `AtMCPoint` objects directly. + +In this tree, `AtSimpleSimulation` uses straight-line propagation. Future versions may extend this step to non-linear tracks and more general transport models. + +### Shared Downstream Stages + +Once `AtMCPoint` objects exist, the remaining stages are shared: + +- `AtClusterizeTask` converts energy deposits into ionization-electron clusters represented as `AtSimulatedPoint` +- `AtPulseTask` drifts those electrons to the pad plane and produces simulated traces in `AtRawEvent` + ## Key Runtime Objects - `AtVertexPropagator` shared simulation-side state between generators and downstream logic - `AtMCTrack` simulated particle tracks -- `AtTpcPoint` - detector-specific transport hits consumed by current digitization tasks - `AtMCPoint` - generic base-family concept for detector points, but usually not the concrete branch a digitization task reads + MC-point type used by digitization logic; also the concrete hit type written by `AtSimpleSimulation` - `AtRawEvent` final simulated raw traces used by reconstruction See [data-model.md](../reference/data-model.md) for the object-level view and [branch-io-contracts.md](../reference/branch-io-contracts.md) for task branch names. -## Running a Simulation +## Event Structure + +The FairRoot/VMC generator path represents each physical beam-induced event as two consecutive FairRoot events: + +- Even-indexed event (0, 2, 4, …): beam phase — the beam particle traverses the detector +- Odd-indexed event (1, 3, 5, …): reaction phase — the reaction products are transported + +Events 0+1 form one complete beam-induced event; events 2+3 form the next, and so on. Code that loops over events or checks event indices must account for this pairing. + +**Track IDs:** Within each FairRoot event the beam particle always has `fTrackID = 0`. Reaction products receive subsequent IDs. + +## Required Pieces -Simulation macros live mostly in `macro/Simulation/` and `macro/examples/`. A typical macro: +For the FairRoot/VMC path, a simulation run needs: -1. Creates a `FairRunSim` -2. Attaches the detector geometry -3. Configures a `FairPrimaryGenerator` with one or more `AtReactionGenerator` subclasses -4. Adds `AtClusterizeTask` and `AtPulseTask` -5. Calls `fRun->Run(nEvents)` +- detector geometry +- a configured `FairPrimaryGenerator` with one or more `AtReactionGenerator` subclasses +- the digitization stages `AtClusterizeTask` and `AtPulseTask` +- the experiment parameter set used by digitization -The output ROOT file contains MC truth alongside detector-specific point branches and the simulated `AtRawEvent`, so it can be fed directly into the reconstruction pipeline. +For the `AtSimpleSimulation` path, the run replaces VMC transport with: -Geometry comes from `geometry/`; runtime parameters typically come from `parameters/`. +- an `AtSimpleSimulation` instance +- one or more configured `AtTools::AtELossModel` instances +- the same downstream digitization stages if `AtRawEvent` output is needed -See [generators.md](generators.md) for generator behavior and [macro-cookbook.md](../reference/macro-cookbook.md) for example macros. +See [generators.md](generators.md) for generator behavior and [energy-loss.md](energy-loss.md) for the model layer used by `AtSimpleSimulation`. diff --git a/docs/subsystems/visualization.md b/docs/subsystems/visualization.md new file mode 100644 index 000000000..5185113e3 --- /dev/null +++ b/docs/subsystems/visualization.md @@ -0,0 +1,112 @@ +# Data Visualization + +ATTPCROOT includes an interactive event display built on ROOT's TEve framework. The current API is `AtViewerManager` (in `AtEventDisplay/`). + +## Architecture + +`AtViewerManager` is a singleton that owns a `FairRunAna` internally and steps through events. It is configured by attaching **tabs** (display panels) and optionally **tasks** (FairTask subclasses for online re-analysis). Branch names for the three standard data streams — `AtRawEvent`, `AtEvent`, and `AtPatternEvent` — are selected at runtime from dropdown menus in the GUI sidebar. + +The viewer requires a pad map (`AtMap`) loaded from an experiment-specific XML mapping file. + +## Basic Setup + +```cpp +// Create the FairRoot run as usual +FairRunAna *fRun = new FairRunAna(); +fRun->SetSource(new FairFileSource(inputFile)); +fRun->SetSink(new FairRootFileSink(outputFile)); +fRun->SetGeomFile(geoFile); + +// Load the pad mapping +TString mapDir = TString(getenv("VMCWORKDIR")) + "/scripts/e12014_pad_mapping.xml"; +auto fMap = std::make_shared(); +fMap->ParseXMLMap(mapDir.Data()); + +// Create viewer and add tabs +AtViewerManager *eveMan = new AtViewerManager(fMap); + +auto tabMain = std::make_unique(); +tabMain->SetMultiHit(100); +eveMan->AddTab(std::move(tabMain)); + +eveMan->Init(); +``` + +## Adding FairTasks for Online Re-analysis + +`AddTask` registers a `FairTask` that runs each time an event is loaded or re-analyzed (e.g. when PSA parameters are changed interactively): + +```cpp +eveMan->AddTask(new AtPSAtask(std::move(psa))); +``` + +Combined with the PSA sidebar addons (`AtSidebarPSA`, `AtSidebarPSADeconv`, etc.), this allows live adjustment of PSA parameters and immediate re-display without restarting the macro. + +## Available Tabs + +### `AtTabMain` +The primary view. Shows a 3D hit display (TEve), a pad plane heat map, and the waveform of the currently selected pad. Selecting a pad in the pad plane updates the waveform display. + +```cpp +auto tab = std::make_unique(); +tab->SetThreshold(20); // minimum charge to draw a hit +tab->SetMultiHit(100); // max hits in a single pad before that pad is suppressed +tab->SetDrawProjection(true); // draw AtPattern projections on the pad plane +tab->SetHitAttributes(TAttMarker(kPink, kFullDotMedium, 1)); +``` + +### `AtTabBraggCurve` +Extends `AtTabMain`. Adds a Bragg curve panel (dE/dx vs. range) alongside the 3D view, populated from `AtTrack::BraggCurve` data. + +### `AtTabPad` +A configurable grid of pad trace panels. Each cell in the grid can independently display a different data source for the currently selected pad. + +```cpp +// 2-row × 4-column grid +auto tab = std::make_unique(2, 4, "MyTraces"); +tab->DrawADC(0, 0); // pedestal-subtracted ADC in cell (0,0) +tab->DrawRawADC(0, 1); // raw ADC +tab->DrawArrayAug("Qreco", 0, 2); // named array augment (e.g. deconvolved charge) +tab->DrawAuxADC("MeshSignal", 0, 3); // auxiliary pad by name +tab->DrawFPN(1, 0); // fixed pattern noise trace +tab->DrawGenTrace(3, 1, 1); // generic trace by ID +tab->DrawHits(0, 0); // overlay hit markers on the ADC trace +eveMan->AddTab(std::move(tab)); +``` + +### `AtTabMacro` +A canvas-based tab driven entirely by user-supplied lambda functions. Callbacks fire at three granularities: + +```cpp +auto tab = std::make_unique(1, 2, "Custom"); + +// Called once when the tree is opened +tab->SetDrawTreeFunction([](AtTabInfo *info) { /* fill histograms */ }, 0, 0); + +// Called each time an event is loaded +tab->SetDrawEventFunction([](AtTabInfo *info) { /* update per-event plots */ }, 0, 1); + +// Called each time the selected pad changes +tab->SetDrawPadFunction([](AtTabInfo *info, Int_t padNum) { /* pad-level plots */ }, 0, 0); + +eveMan->AddTab(std::move(tab)); +``` + +### `AtTabFission` +Extends `AtTabMain`. Adds hit overlays for fission fragment analysis (uncorrected and corrected hit sets). Exposes `GetFissionBranch()` for branch selection. + +### `AtTabEnergyLoss` +Canvas-based tab that calculates dE/dx for two fission fragments from an `AtFissionEvent`. Stacked histograms, charge-sum and Gaussian-fit methods, and the ratio between the two fragments. + +```cpp +AtTabFission *fissionTab = ...; // add fission tab first +auto elTab = std::make_unique(fissionTab->GetFissionBranch()); +eveMan->AddTab(std::move(elTab)); +``` + +## Notes + +- **Branch selection is in the GUI**, not in the macro. The sidebar dropdowns populate from whatever branches are present in the input file. +- `AtViewerManager::Instance()` returns the singleton after construction. +- Requires a display (X11 or native windowing); cannot run headlessly. +- Pad mapping XML files live in `scripts/`. The correct file is experiment-specific. diff --git a/docs/tooling/testing.md b/docs/tooling/testing.md index be5a8a78e..fefe3264d 100644 --- a/docs/tooling/testing.md +++ b/docs/tooling/testing.md @@ -47,6 +47,27 @@ These are listed and run separately from unit tests. Use them for: - simulation plus reconstruction flows - experiment-facing regression scenarios +### Prerequisites + +Integration tests require geometry ROOT files that are **not** committed to the repo. Generate them first: + +```bash +macro/tests/generateGeometry.sh +``` + +Some tests also require external data files that are too large for the repo: + +- **GADGETII tests**: require the fishtank dataset (available on FRIB/NSCL cluster only) +- **SpecMAT tests**: require `TTreesGETrun_9993.root` placed in `macro/tests/SpecMAT/data/` + +### Running All Integration Tests + +```bash +macro/tests/runAllTest.sh +``` + +Results are logged in each test subdirectory as `test.log`. + See [macro-cookbook.md](macro-cookbook.md) for curated macro starting points. ## Continuous Integration From 055b256d58e378229c42da0ff338182052c6f346 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Wed, 1 Apr 2026 20:20:17 -0400 Subject: [PATCH 08/10] Fix errors in doc --- docs/contributing/code-style.md | 8 +++--- docs/contributing/guide.md | 1 + docs/reference/parameters.md | 46 ++++++++++++++++---------------- docs/subsystems/visualization.md | 4 +-- 4 files changed, 31 insertions(+), 28 deletions(-) diff --git a/docs/contributing/code-style.md b/docs/contributing/code-style.md index e0062ee0b..f62af2648 100644 --- a/docs/contributing/code-style.md +++ b/docs/contributing/code-style.md @@ -25,14 +25,16 @@ Follow ROOT naming conventions: ### ClassDef / ClassImp -Use `ClassDef` and `ClassImp` **only** for `TObject` subclasses. When the memory layout changes, increment the version number in `ClassDef`: +Use `ClassDef` and `ClassImp` only for classes that need ROOT-generated class metadata or streamer support. + +For reflection-only classes, prefer a `-!` LinkDef entry. Use `+` only for classes written to disk. + +When a streamed class changes memory layout, increment the version number in `ClassDef`: ```cpp ClassDef(AtMyData, 2); // bump when adding/removing/reordering members ``` -Non-persisted classes (tasks, algorithms, models) do not need `ClassDef`. - ## Formatting The project uses `clang-format-17` with the repo-root `.clang-format`. diff --git a/docs/contributing/guide.md b/docs/contributing/guide.md index 73ef7a25a..7713eba04 100644 --- a/docs/contributing/guide.md +++ b/docs/contributing/guide.md @@ -22,6 +22,7 @@ Double_t fValue; ## ROOT Dictionaries Every class with `ClassDef` must appear in the appropriate `*LinkDef.h`. +Classes listed in `*LinkDef.h` do not automatically need `ClassDef` or `ClassImp`. Preferred suffix rules for new or edited code: diff --git a/docs/reference/parameters.md b/docs/reference/parameters.md index cacd18d88..326db7f70 100644 --- a/docs/reference/parameters.md +++ b/docs/reference/parameters.md @@ -12,10 +12,12 @@ rtdb->setFirstInput(parIo); rtdb->getContainer("AtDigiPar"); ``` -Parameter files live in `parameters/`. Each file typically contains an `[AtDigiPar]` section followed by key-value pairs. Pick the file that matches the experiment or create a new one by copying an existing one. +Parameter files live in `parameters/`. Each file typically contains an `[AtDigiPar]` section followed by key-value pairs. Pick the file that matches the experiment or create a new one by copying a current experiment file such as `parameters/ATTPC.e20020.par`. ## AtDigiPar Parameters +`AtDigiPar::getParams()` requires the following keys: + ### Electromagnetic Fields | Parameter | Type | Unit | Description | @@ -29,7 +31,6 @@ Parameter files live in `parameters/`. Each file typically contains an `[AtDigiP |-----------|------|------|-------------| | `TBEntrance` | Int_t | time buckets | Beam entrance position at detector entrance | | `ZPadPlane` | Double_t | mm | Position of the micromegas pad plane (detector length) | -| `DriftLength` | Double_t | mm | Full drift length (AT-TPC: 1000, prototype: 500) | ### Gas Physics @@ -38,8 +39,8 @@ Parameter files live in `parameters/`. Each file typically contains an `[AtDigiP | `DriftVelocity` | Double_t | cm/µs | Electron drift velocity | | `EIonize` | Double_t | eV | Effective ionization energy of the fill gas | | `Fano` | Double_t | — | Fano factor of the fill gas | -| `CoefL` | Double_t | cm^−½ | Longitudinal diffusion coefficient | -| `CoefT` | Double_t | cm^−½ | Transverse diffusion coefficient | +| `CoefL` | Double_t | cm²/µs | Longitudinal diffusion coefficient | +| `CoefT` | Double_t | cm²/µs | Transverse diffusion coefficient | | `GasPressure` | Double_t | torr | Gas pressure | | `Density` | Double_t | kg/m³ | Gas density | @@ -51,7 +52,6 @@ Parameter files live in `parameters/`. Each file typically contains an `[AtDigiP | `GETGain` | Double_t | fC | Gain from GET electronics | | `PeakingTime` | Int_t | ns | Electronic response peaking time | | `Gain` | Double_t | — | Average micromegas amplification factor | -| `NumTbs` | Int_t | — | Number of time buckets | ## Accessing Parameters in Code @@ -69,25 +69,25 @@ double driftVel = fPar->GetDriftVelocity(); // cm/µs double bField = fPar->GetBField(); // T ``` -## Example Parameter File Entry +## Example `AtDigiPar` Block -``` +```text [AtDigiPar] -EField:Double_t 5000 # V/m -BField:Double_t 2 # Tesla -TBEntrance:Int_t 280 # time bucket at detector entrance -ZPadPlane:Double_t 1000 # mm -DriftVelocity:Double_t 5.00 # cm/us -EIonize:Double_t 15.603 # eV (He+CO2 example) -Fano:Double_t 0.22 -CoefL:Double_t 0.025 # cm^-0.5 -CoefT:Double_t 0.010 # cm^-0.5 -GasPressure:Double_t 100 # torr -Density:Double_t 0.0738 # kg/m3 -SamplingRate:Int_t 12.5 # MHz -Gain:Double_t 100.0 -NumTbs:Int_t 512 -PeakingTime:Int_t 117 # ns +BField:Double_t 2.0 +EField:Double_t 5000 +TBEntrance:Int_t 445 +ZPadPlane:Double_t 1000.0 +EIonize:Double_t 15.603 +Fano:Double_t 0.24 +CoefL:Double_t 0.0003 +CoefT:Double_t 0.0003 +GasPressure:Double_t 600 +Density:Double_t 0.197 +DriftVelocity:Double_t 0.87 +Gain:Double_t 1000 +SamplingRate:Int_t 3 +GETGain:Double_t 1000 +PeakingTime:Int_t 720 ``` -See `parameters/AT.parameters.par` for a complete reference example and `parameters/` for all experiment-specific files. +Copy a complete file from `parameters/` and edit the values for the target experiment. diff --git a/docs/subsystems/visualization.md b/docs/subsystems/visualization.md index 5185113e3..529773ee4 100644 --- a/docs/subsystems/visualization.md +++ b/docs/subsystems/visualization.md @@ -4,14 +4,14 @@ ATTPCROOT includes an interactive event display built on ROOT's TEve framework. ## Architecture -`AtViewerManager` is a singleton that owns a `FairRunAna` internally and steps through events. It is configured by attaching **tabs** (display panels) and optionally **tasks** (FairTask subclasses for online re-analysis). Branch names for the three standard data streams — `AtRawEvent`, `AtEvent`, and `AtPatternEvent` — are selected at runtime from dropdown menus in the GUI sidebar. +`AtViewerManager` is a singleton GUI/controller. It uses `FairRunAna::Instance()` in `AddTask()` and `Init()`, so the macro must create and configure `FairRunAna` before constructing the viewer. Attach **tabs** (display panels) and optionally **tasks** (FairTask subclasses for online re-analysis). Branch names for the three standard data streams — `AtRawEvent`, `AtEvent`, and `AtPatternEvent` — are selected at runtime from dropdown menus in the GUI sidebar. The viewer requires a pad map (`AtMap`) loaded from an experiment-specific XML mapping file. ## Basic Setup ```cpp -// Create the FairRoot run as usual +// Create and configure the FairRoot run before constructing the viewer FairRunAna *fRun = new FairRunAna(); fRun->SetSource(new FairFileSource(inputFile)); fRun->SetSink(new FairRootFileSink(outputFile)); From 70036e13ea1b6e6481c1d2ce6d5a415cf393d392 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Wed, 1 Apr 2026 21:47:39 -0400 Subject: [PATCH 09/10] Remove unused files --- ...is-implementing-a-kalman-filter-for-us.txt | 867 ------------------ .claude/plans/keen-marinating-wind.md | 295 ------ .claude/resume-ukf-task-integration.md | 38 - 3 files changed, 1200 deletions(-) delete mode 100644 .claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt delete mode 100644 .claude/plans/keen-marinating-wind.md delete mode 100644 .claude/resume-ukf-task-integration.md diff --git a/.claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt b/.claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt deleted file mode 100644 index cf1499347..000000000 --- a/.claude/2026-03-03-224207-this-branch-is-implementing-a-kalman-filter-for-us.txt +++ /dev/null @@ -1,867 +0,0 @@ - -╭─── Claude Code v2.1.63 ──────────────────────────────────────────────────────╮ -│ │ Tips for getting │ -│ Welcome back Adam! │ started │ -│ │ Run /init to create a … │ -│ ▐▛███▜▌ │ ─────────────────────── │ -│ ▝▜█████▛▘ │ Recent activity │ -│ ▘▘ ▝▝ │ No recent activity │ -│ │ │ -│ Sonnet 4.6 · Claude Pro · aanthony3@gmail.com's │ │ -│ Organization │ │ -│ ~/fair_install/ATTPCROOTv2 │ │ -╰──────────────────────────────────────────────────────────────────────────────╯ - -❯ This branch is implementing a Kalman Filter for use in TPCs. You can see - tests showing how it is being used to fit tracks in @macro/tests/UKF/. - The class itself is located in - @AtReconstruction/AtFitter/OpenKF/kalman_filter/TrackFitterUKF* and it - depends on the propagator defined in @AtTools/AtPropagator*. This fitter - needs to be integrated into the task system of the codebase. The output - will be an AtFittedTrack which has a proposed class in the PR #254 which - is to merge the branch RealAurio:fitterRefactor. Please make a plan to - finalize this refactor and implement a task that will run the kalman - filter on tracks. The AtFitter task should share a base with the other - fitter tasks. At least the GenFit version, but also the MCFitter if it - combines well. - ⎿  Listed directory macro/tests/UKF/ - -● 3 Explore agents finished (ctrl+o to expand) - ├─ Explore UKF fitter and propagator classes · 12 tool uses · 55.6k tokens - │ ⎿  Done - ├─ Explore existing fitter tasks and task system · 25 tool uses · 50.2k - tokens - │ ⎿  Done - └─ Explore AtFittedTrack and PR #254 branch · 43 tool uses · 55.1k tokens - ⎿  Done - -● Let me read the critical files directly before designing the plan. - -● Read 6 files (ctrl+o to expand) - -● Now let me look at the CMakeLists files and the existing OpenKFTest to - understand build system integration: - -● Searched for 2 patterns, read 1 file (ctrl+o to expand) - -● Let me check the AtHitCluster and AtTrack data to understand what information - is available for UKF: - -● Searched for 2 patterns, read 1 file (ctrl+o to expand) - -● Now let me look at TrackFitterUKF's physics-specific interface in detail: - -● Searched for 1 pattern, read 1 file (ctrl+o to expand) - -● Now let me quickly check the ELoss model setup and check what the AtGenfit - does for momentum estimation: - -● Searched for 4 patterns (ctrl+o to expand) - -● Before finalizing the plan, I need to clarify a few key decisions: - ⎿  Tool use rejected with user message: The user wants to clarify these - questions. - This means they may have additional information, context or questions - for you. - Take their response into account and then reformulate the questions if - appropriate. - Start by asking them what they would like to clarify. - - Questions asked: - - "The fitterRefactor PR (#254) changes ~100+ files (namespace rename - AtFITTER→EventFit, AtFitter interface change from ProcessTracks() to - FitEvent()/GetFittedTrack(), new AtFittedTrack nested structs, removes PRA - data). How much of this refactor should we adopt?" - (No answer provided) - - "The UKF needs an initial momentum magnitude as a seed. What is the - primary strategy for providing this?" - (No answer provided) - - "What UKF output should be stored in AtFittedTrack beyond the standard - vertex state (energy, theta, phi, position)?" - (No answer provided) - -● What would you like to clarify about the questions I asked? - -✻ Baked for 4m 23s - -❯ I think you are comparing the wrong things. main is not the main branch, - develop is. Main is very old and should be ignored. That PR only touches - 36 files, not 100. Those changes are all important to understand to - construct a plan against them. To set the initial momentum, we will use a - circle fit like we do in the getfit case (PRA). Give me your options for - the last question again. - -● Good point - let me look at the fitterRefactor changes properly (against - develop) and how GenFit uses the circle fit for initial momentum. - -● Bash(git log --oneline RealAurio/fitterRefactor ^develop | head -20) - ⎿  11b7a6f9 Removed stuff from AtFitter that ended up being too case spec - ific to make sense to keep in the base class. Also fixed some serious - bugs that I found while actually implementing a fitter class in the br - … +16 lines (ctrl+o to expand) - -● Bash(git diff --name-only RealAurio/fitterRefactor develop -- 2>/dev/null | - head -50) - ⎿  AtData/AtDataLinkDef.h - AtData/AtFitMetadata.cxx - AtData/AtFitMetadata.h - … +27 lines (ctrl+o to expand) - -● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtFitter.h) - ⎿  #ifndef AtFITTER_H - #define AtFITTER_H - … +40 lines (ctrl+o to expand) - -● Bash(git show RealAurio/fitterRefactor:AtData/AtFittedTrack.h) - ⎿  #ifndef ATFITTEDTRACK_H - #define ATFITTEDTRACK_H - … +162 lines (ctrl+o to expand) - -● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitterTask.h && echo - "---" && git show - RealAurio/fitterRefactor:AtReconstruction/AtFitterTask.cxx) - ⎿  /********************************************************************* - * Fitter Task AtFitterTask.hh * - * Author: Y. Ayyad ayyadlim@frib.msu.edu * - … +95 lines (ctrl+o to expand) - -● Bash(git show RealAurio/fitterRefactor:AtData/AtFitTrackMetadata.h && echo - "---" && git show RealAurio/fitterRefactor:AtData/AtFitMetadata.h) - ⎿  #ifndef ATFITTRACKMETADATA_H - #define ATFITTRACKMETADATA_H - … +39 lines (ctrl+o to expand) - -● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtFitter.cxx) - ⎿  #include "AtFitter.h" - - #include "AtPatternEvent.h" - … +25 lines (ctrl+o to expand) - -● Searched for 2 patterns, read 1 file (ctrl+o to expand) - -● I now have the full picture. Here's the re-ask on output storage: - -● User answered Claude's questions: - ⎿  · What UKF output should be stored in AtFittedTrack beyond the - standard vertex state (Kinematics + Vertex + TrackProperties already - in the refactored class)? → Vertex + smoothed positions (Recommended) - -● Bash(git show RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtGenfit.h | - head -60) - ⎿  #ifndef ATGENFIT_H - #define ATGENFIT_H - … +36 lines (ctrl+o to expand) - -● Bash(git show RealAurio/fitterRefactor:AtData/AtTrack.h | grep -n - "smoothed\|Smooth\|UKF\|fSmooth\|fFitted" | head -20) - ⎿  (No output) - -● Bash(git show - RealAurio/fitterRefactor:AtReconstruction/AtFitter/AtFitterOld.h) - ⎿  #ifndef AtFITTEROLD_H - #define AtFITTEROLD_H - … +37 lines (ctrl+o to expand) - -● Bash(git show RealAurio/fitterRefactor:AtData/AtTrackingEventOld.h | head -30) - ⎿  #ifndef AtTRACKINGEVENTOLD_H - #define AtTRACKINGEVENTOLD_H - … +27 lines (ctrl+o to expand) - -● Searched for 1 pattern (ctrl+o to expand) - -● Now I have everything I need. Let me write the plan: - -● Updated plan - ⎿  /plan to preview - ⎿  User rejected Claude's plan: - ╭─────────────────────────────────────────────────────────────────────────╮ - │ Plan: Integrate UKF Track Fitter into Task System │ - │ │ - │ Context │ - │ │ - │ The UKF track fitter (kf::TrackFitterUKF) is fully implemented as a │ - │ standalone physics class and tested via macros in macro/tests/UKF/. It │ - │ needs to be integrated into the FairTask-based reconstruction pipeline │ - │ so it can be run event-by-event on data. │ - │ │ - │ In parallel, PR #254 (RealAurio/fitterRefactor) proposes a clean │ - │ refactored fitter architecture with a new EventFit::AtFitter base class │ - │ and a redesigned AtFittedTrack. The UKF task should be implemented │ - │ against this new interface rather than the deprecated │ - │ AtFITTER::AtFitter. │ - │ │ - │ Outcome: A complete EventFit::AtFitterUKF class that reads an │ - │ AtPatternEvent (via AtFitterTask), fits each track with the UKF using a │ - │ Brho-based momentum seed, and writes an AtTrackingEvent with │ - │ AtFittedTrack objects containing vertex kinematics and smoothed │ - │ trajectory positions. │ - │ │ - │ --- │ - │ Step 1 — Merge fitterRefactor (PR #254) into OpenKF-Claude │ - │ │ - │ Merge RealAurio/fitterRefactor into the OpenKF-Claude branch. This │ - │ brings in all 36 changed files. Key files introduced/changed: │ - │ │ - │ File: AtData/AtFittedTrack.h/cxx │ - │ Change: New nested-struct design (Kinematics, ParticleInfo, │ - │ TrackProperties) │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtFitMetadata.h/cxx │ - │ Change: New: event-level fit metadata │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtFitTrackMetadata.h/cxx │ - │ Change: New: per-track fit statistics │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtFittedTrackOld.h/cxx │ - │ Change: Backward-compat alias for old AtFittedTrack │ - │ ──────────────────────────────────────── │ - │ File: AtData/CMakeLists.txt │ - │ Change: Adds new data classes │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtDataLinkDef.h │ - │ Change: Adds ROOT dictionary entries │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtFitter.h/cxx │ - │ Change: New EventFit::AtFitter with FitEvent()/GetFittedTrack() │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtFitterOld.h/cxx │ - │ Change: Old AtFITTER::AtFitter kept for backward compat │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtGenfit.h │ - │ Change: Now inherits AtFitterOld (deprecated, not changed) │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtMCFitter.cxx │ - │ Change: Minor updates │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitterTask.h/cxx │ - │ Change: Updated to use EventFit::AtFitter; adds optional │ - │ RawEvent/Event/FitMetadata branches │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtReconstructionLinkDef.h │ - │ Change: Updated LinkDef │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtPatternModification/* │ - │ Change: New BraggCurveFinder task (independent) │ - │ ──────────────────────────────────────── │ - │ File: AtEventDisplay/AtTabs/* │ - │ Change: Tab updates (independent) │ - │ │ - │ Resolve any merge conflicts; most conflict risk is in CMakeLists and │ - │ LinkDef files. │ - │ │ - │ --- │ - │ Step 2 — Add Smoothed Positions to AtFittedTrack::TrackProperties │ - │ │ - │ File: AtData/AtFittedTrack.h │ - │ │ - │ Add to TrackProperties struct: │ - │ std::vector fSmoothedPositions; // UKF-fitted │ - │ positions at each cluster │ - │ │ - │ Add corresponding setter/getter to AtFittedTrack: │ - │ void SetSmoothedPositions(std::vector pos) { │ - │ fTrackProperties.fSmoothedPositions = std::move(pos); } │ - │ const std::vector& GetSmoothedPositions() const { │ - │ return fTrackProperties.fSmoothedPositions; } │ - │ │ - │ Bump ClassDef(AtFittedTrack, 3). │ - │ │ - │ --- │ - │ Step 3 — Implement AtFitterUKF │ - │ │ - │ New files: │ - │ - AtReconstruction/AtFitter/AtFitterUKF.h │ - │ - AtReconstruction/AtFitter/AtFitterUKF.cxx │ - │ │ - │ Header │ - │ │ - │ namespace EventFit { │ - │ │ - │ class AtFitterUKF : public AtFitter { │ - │ public: │ - │ AtFitterUKF(double charge, double mass_MeV, │ - │ std::unique_ptr elossModel); │ - │ │ - │ void SetBField(ROOT::Math::XYZVector bField); // Tesla │ - │ void SetEField(ROOT::Math::XYZVector eField); // V/m (default 0) │ - │ void SetUKFParameters(double alpha, double beta, double kappa); │ - │ void SetMeasurementSigma(double sigma_mm); // default 1.0 mm │ - │ void SetMomentumSigmaFrac(double frac); // fractional, default │ - │ 0.1 │ - │ void SetMinClusters(int n); // default 3 │ - │ │ - │ void Init() override {} │ - │ │ - │ protected: │ - │ AtFittedTrack *GetFittedTrack(AtTrack *track, │ - │ AtFitMetadata *fitMetadata = nullptr, │ - │ AtRawEvent *rawEvent = nullptr, │ - │ AtEvent *event = nullptr) override; │ - │ private: │ - │ // Helpers │ - │ std::unique_ptr CreateUKF() const; │ - │ ROOT::Math::XYZPoint GetInitialPosition(AtTrack *track) const; │ - │ ROOT::Math::XYZVector GetInitialMomentum(AtTrack *track) const; │ - │ TMatrixD GetInitialCovariance(double p_mag_MeV) const; │ - │ TMatrixD GetMeasCovariance() const; │ - │ double GetBrho(AtTrack *track) const; │ - │ │ - │ // Particle / detector config │ - │ double fCharge; // charge in │ - │ units of e │ - │ double fMass_MeV; // rest mass │ - │ in MeV/c² │ - │ ROOT::Math::XYZVector fBField{0, 0, 2.85}; // Tesla │ - │ ROOT::Math::XYZVector fEField{0, 0, 0}; // V/m │ - │ std::unique_ptr fELossModel; │ - │ │ - │ // UKF tuning │ - │ double fAlpha{1e-3}, fBeta{2.0}, fKappa{0.0}; │ - │ double fMeasSigma_mm{1.0}; │ - │ double fMomSigmaFrac{0.1}; │ - │ int fMinClusters{3}; │ - │ │ - │ AtTools::AtKinematics fKinematics; // for GetMomFromBrho │ - │ }; │ - │ │ - │ } // namespace EventFit │ - │ │ - │ Implementation logic in GetFittedTrack() │ - │ │ - │ 1. Get clusters: auto *clusters = track->GetHitClusterArray(); │ - │ track->SortClusterHitArrayZ(); │ - │ Skip track if clusters->size() < fMinClusters. │ - │ 2. Seed momentum from circle fit (same as GenFit): │ - │ radius = track->GetGeoRadius() / 1000.0 // mm → m │ - │ theta = track->GetGeoTheta() │ - │ phi = track->GetGeoPhi() │ - │ brho = fBField.Z() * radius / sin(theta) // T·m │ - │ [p_GeV, KE_GeV] = fKinematics.GetMomFromBrho(mass_amu, charge, brho) │ - │ p_MeV = p_GeV * 1000.0 │ - │ 3. Initial position: clusters->front().GetPosition() (in mm). │ - │ 4. Initial momentum direction: unit vector from (theta, phi); magnitude │ - │ = p_MeV. │ - │ 5. Create and configure UKF: │ - │ auto ukf = CreateUKF(); // new kf::TrackFitterUKF with │ - │ AtPropagator+AtRK4AdaptiveStepper │ - │ ukf->setParameters(fAlpha, fBeta, fKappa); │ - │ ukf->SetInitialState(initialPos, initialMom, │ - │ GetInitialCovariance(p_MeV)); │ - │ ukf->SetMeasCov(GetMeasCovariance()); │ - │ ukf->fEnableEnStraggling = true; │ - │ 6. Filter loop: │ - │ for (size_t i = 1; i < clusters->size(); ++i) { │ - │ XYZPoint measPt = clusters->at(i).GetPosition(); // mm │ - │ ukf->predictUKF(measPt); │ - │ ukf->correctUKF(measPt); │ - │ } │ - │ ukf->smoothUKF(); │ - │ 6. Wrap in try/catch for Cholesky failures → mark fitConverged = false. │ - │ 7. Build AtFittedTrack: │ - │ - From smoothedStates[0] = [x, y, z, p, theta_s, phi_s]: │ - │ - KE = sqrt(p*p + m*m) - m │ - │ - SetKinematics(KE, theta_s, phi_s) │ - │ - SetVertex(XYZVector{x, y, z}) │ - │ - SetParticleInfo(pdg_string, charge, mass_amu) │ - │ - Extract smoothed positions from all states → │ - │ SetSmoothedPositions(...) │ - │ - Compute trackLength as sum of distances between consecutive │ - │ smoothed positions │ - │ - SetTrackPropertiesStruct(...) with initialPosition, trackLength, │ - │ trackPoints = clusters->size() │ - │ - Create AtFitTrackMetadata: set trackID, fitConverged; chi² = sum of │ - │ squared smoothed residuals / ndf; SetTrackMetadata() │ - │ 8. Optionally store metadata in AtFitMetadata (if non-null) via │ - │ SetTrackMetadatasVector(trackID, ...). │ - │ │ - │ --- │ - │ Step 4 — Update Build System │ - │ │ - │ File: AtReconstruction/CMakeLists.txt │ - │ - Add AtFitter/AtFitterUKF.cxx to SRCS │ - │ │ - │ File: AtReconstruction/AtReconstructionLinkDef.h │ - │ - Add #pragma link C++ class EventFit::AtFitterUKF+; │ - │ │ - │ --- │ - │ Step 5 — Write Usage Macro │ - │ │ - │ New file: macro/tests/UKF/UKFTask.C │ - │ │ - │ Demonstrates: │ - │ auto eloss = std::make_unique(3.553e-5); │ - │ eloss->SetProjectile(1, 1, 1); │ - │ eloss->SetMaterial({{1, 1, 1}}); │ - │ │ - │ auto fitter = std::make_unique(1.0, 938.272, │ - │ std::move(eloss)); │ - │ fitter->SetBField({0, 0, 2.85}); │ - │ fitter->SetUKFParameters(1e-3, 2, 0); │ - │ fitter->SetMeasurementSigma(1.0); │ - │ │ - │ auto task = new AtFitterTask(std::move(fitter)); │ - │ task->SetPersistence(kTRUE); │ - │ run->AddTask(task); │ - │ │ - │ --- │ - │ MCFitter Assessment │ - │ │ - │ AtMCFitter uses stochastic multi-iteration Monte Carlo optimization and │ - │ produces AtMCResult arrays — not per-track AtFittedTrack objects. Its │ - │ Exec(AtPatternEvent&) pattern is fundamentally different from │ - │ GetFittedTrack(AtTrack*). It does not combine well with the new │ - │ EventFit::AtFitter base and should remain separate (AtMCFitterTask). │ - │ │ - │ --- │ - │ Verification │ - │ │ - │ 1. Build: cmake --build build -j$(nproc) — must compile without errors. │ - │ 2. Existing macro: Run macro/tests/UKF/UKFSingleTrack.C to confirm │ - │ TrackFitterUKF still works standalone. │ - │ 3. New macro: Run macro/tests/UKF/UKFTask.C end-to-end on a simulation │ - │ file; inspect the AtTrackingEvent branch in the output ROOT file — │ - │ verify AtFittedTrack contains non-trivial kinematics and smoothed │ - │ positions. │ - │ 4. Check merge: Verify AtGenfit still compiles (it uses AtFitterOld │ - │ which is unchanged). │ - ╰─────────────────────────────────────────────────────────────────────────╯ - -❯ To this plan, I think we should emphasize the need for proper tests using - the gtest stuff we use other places in the code-base. - -● Explore(Explore gtest structure in codebase) - ⎿  Done (30 tool uses · 52.4k tokens · 1m 4s) - (ctrl+o to expand) - -● Good. Let me also read the existing UKF test file to understand what's already - covered: - -● Read 1 file (ctrl+o to expand) - -● Now I have a complete picture of the test infrastructure. Let me update the - plan with a proper testing section: - -● Updated plan - ⎿  /plan to preview - ⎿  User rejected Claude's plan: - ╭─────────────────────────────────────────────────────────────────────────╮ - │ Plan: Integrate UKF Track Fitter into Task System │ - │ │ - │ Context │ - │ │ - │ The UKF track fitter (kf::TrackFitterUKF) is fully implemented as a │ - │ standalone physics class and tested via macros in macro/tests/UKF/. It │ - │ needs to be integrated into the FairTask-based reconstruction pipeline │ - │ so it can be run event-by-event on data. │ - │ │ - │ In parallel, PR #254 (RealAurio/fitterRefactor) proposes a clean │ - │ refactored fitter architecture with a new EventFit::AtFitter base class │ - │ and a redesigned AtFittedTrack. The UKF task should be implemented │ - │ against this new interface rather than the deprecated │ - │ AtFITTER::AtFitter. │ - │ │ - │ Outcome: A complete EventFit::AtFitterUKF class that reads an │ - │ AtPatternEvent (via AtFitterTask), fits each track with the UKF using a │ - │ Brho-based momentum seed, and writes an AtTrackingEvent with │ - │ AtFittedTrack objects containing vertex kinematics and smoothed │ - │ trajectory positions. │ - │ │ - │ --- │ - │ Step 1 — Merge fitterRefactor (PR #254) into OpenKF-Claude │ - │ │ - │ Merge RealAurio/fitterRefactor into the OpenKF-Claude branch. This │ - │ brings in all 36 changed files. Key files introduced/changed: │ - │ │ - │ File: AtData/AtFittedTrack.h/cxx │ - │ Change: New nested-struct design (Kinematics, ParticleInfo, │ - │ TrackProperties) │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtFitMetadata.h/cxx │ - │ Change: New: event-level fit metadata │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtFitTrackMetadata.h/cxx │ - │ Change: New: per-track fit statistics │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtFittedTrackOld.h/cxx │ - │ Change: Backward-compat alias for old AtFittedTrack │ - │ ──────────────────────────────────────── │ - │ File: AtData/CMakeLists.txt │ - │ Change: Adds new data classes │ - │ ──────────────────────────────────────── │ - │ File: AtData/AtDataLinkDef.h │ - │ Change: Adds ROOT dictionary entries │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtFitter.h/cxx │ - │ Change: New EventFit::AtFitter with FitEvent()/GetFittedTrack() │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtFitterOld.h/cxx │ - │ Change: Old AtFITTER::AtFitter kept for backward compat │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtGenfit.h │ - │ Change: Now inherits AtFitterOld (deprecated, not changed) │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitter/AtMCFitter.cxx │ - │ Change: Minor updates │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtFitterTask.h/cxx │ - │ Change: Updated to use EventFit::AtFitter; adds optional │ - │ RawEvent/Event/FitMetadata branches │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtReconstructionLinkDef.h │ - │ Change: Updated LinkDef │ - │ ──────────────────────────────────────── │ - │ File: AtReconstruction/AtPatternModification/* │ - │ Change: New BraggCurveFinder task (independent) │ - │ ──────────────────────────────────────── │ - │ File: AtEventDisplay/AtTabs/* │ - │ Change: Tab updates (independent) │ - │ │ - │ Resolve any merge conflicts; most conflict risk is in CMakeLists and │ - │ LinkDef files. │ - │ │ - │ --- │ - │ Step 2 — Add Smoothed Positions to AtFittedTrack::TrackProperties │ - │ │ - │ File: AtData/AtFittedTrack.h │ - │ │ - │ Add to TrackProperties struct: │ - │ std::vector fSmoothedPositions; // UKF-fitted │ - │ positions at each cluster │ - │ │ - │ Add corresponding setter/getter to AtFittedTrack: │ - │ void SetSmoothedPositions(std::vector pos) { │ - │ fTrackProperties.fSmoothedPositions = std::move(pos); } │ - │ const std::vector& GetSmoothedPositions() const { │ - │ return fTrackProperties.fSmoothedPositions; } │ - │ │ - │ Bump ClassDef(AtFittedTrack, 3). │ - │ │ - │ --- │ - │ Step 3 — Implement AtFitterUKF │ - │ │ - │ New files: │ - │ - AtReconstruction/AtFitter/AtFitterUKF.h │ - │ - AtReconstruction/AtFitter/AtFitterUKF.cxx │ - │ │ - │ Header │ - │ │ - │ namespace EventFit { │ - │ │ - │ class AtFitterUKF : public AtFitter { │ - │ public: │ - │ AtFitterUKF(double charge, double mass_MeV, │ - │ std::unique_ptr elossModel); │ - │ │ - │ void SetBField(ROOT::Math::XYZVector bField); // Tesla │ - │ void SetEField(ROOT::Math::XYZVector eField); // V/m (default 0) │ - │ void SetUKFParameters(double alpha, double beta, double kappa); │ - │ void SetMeasurementSigma(double sigma_mm); // default 1.0 mm │ - │ void SetMomentumSigmaFrac(double frac); // fractional, default │ - │ 0.1 │ - │ void SetMinClusters(int n); // default 3 │ - │ │ - │ void Init() override {} │ - │ │ - │ protected: │ - │ AtFittedTrack *GetFittedTrack(AtTrack *track, │ - │ AtFitMetadata *fitMetadata = nullptr, │ - │ AtRawEvent *rawEvent = nullptr, │ - │ AtEvent *event = nullptr) override; │ - │ private: │ - │ // Helpers │ - │ std::unique_ptr CreateUKF() const; │ - │ ROOT::Math::XYZPoint GetInitialPosition(AtTrack *track) const; │ - │ ROOT::Math::XYZVector GetInitialMomentum(AtTrack *track) const; │ - │ TMatrixD GetInitialCovariance(double p_mag_MeV) const; │ - │ TMatrixD GetMeasCovariance() const; │ - │ double GetBrho(AtTrack *track) const; │ - │ │ - │ // Particle / detector config │ - │ double fCharge; // charge in │ - │ units of e │ - │ double fMass_MeV; // rest mass │ - │ in MeV/c² │ - │ ROOT::Math::XYZVector fBField{0, 0, 2.85}; // Tesla │ - │ ROOT::Math::XYZVector fEField{0, 0, 0}; // V/m │ - │ std::unique_ptr fELossModel; │ - │ │ - │ // UKF tuning │ - │ double fAlpha{1e-3}, fBeta{2.0}, fKappa{0.0}; │ - │ double fMeasSigma_mm{1.0}; │ - │ double fMomSigmaFrac{0.1}; │ - │ int fMinClusters{3}; │ - │ │ - │ AtTools::AtKinematics fKinematics; // for GetMomFromBrho │ - │ }; │ - │ │ - │ } // namespace EventFit │ - │ │ - │ Implementation logic in GetFittedTrack() │ - │ │ - │ 1. Get clusters: auto *clusters = track->GetHitClusterArray(); │ - │ track->SortClusterHitArrayZ(); │ - │ Skip track if clusters->size() < fMinClusters. │ - │ 2. Seed momentum from circle fit (same as GenFit): │ - │ radius = track->GetGeoRadius() / 1000.0 // mm → m │ - │ theta = track->GetGeoTheta() │ - │ phi = track->GetGeoPhi() │ - │ brho = fBField.Z() * radius / sin(theta) // T·m │ - │ [p_GeV, KE_GeV] = fKinematics.GetMomFromBrho(mass_amu, charge, brho) │ - │ p_MeV = p_GeV * 1000.0 │ - │ 3. Initial position: clusters->front().GetPosition() (in mm). │ - │ 4. Initial momentum direction: unit vector from (theta, phi); magnitude │ - │ = p_MeV. │ - │ 5. Create and configure UKF: │ - │ auto ukf = CreateUKF(); // new kf::TrackFitterUKF with │ - │ AtPropagator+AtRK4AdaptiveStepper │ - │ ukf->setParameters(fAlpha, fBeta, fKappa); │ - │ ukf->SetInitialState(initialPos, initialMom, │ - │ GetInitialCovariance(p_MeV)); │ - │ ukf->SetMeasCov(GetMeasCovariance()); │ - │ ukf->fEnableEnStraggling = true; │ - │ 6. Filter loop: │ - │ for (size_t i = 1; i < clusters->size(); ++i) { │ - │ XYZPoint measPt = clusters->at(i).GetPosition(); // mm │ - │ ukf->predictUKF(measPt); │ - │ ukf->correctUKF(measPt); │ - │ } │ - │ ukf->smoothUKF(); │ - │ 6. Wrap in try/catch for Cholesky failures → mark fitConverged = false. │ - │ 7. Build AtFittedTrack: │ - │ - From smoothedStates[0] = [x, y, z, p, theta_s, phi_s]: │ - │ - KE = sqrt(p*p + m*m) - m │ - │ - SetKinematics(KE, theta_s, phi_s) │ - │ - SetVertex(XYZVector{x, y, z}) │ - │ - SetParticleInfo(pdg_string, charge, mass_amu) │ - │ - Extract smoothed positions from all states → │ - │ SetSmoothedPositions(...) │ - │ - Compute trackLength as sum of distances between consecutive │ - │ smoothed positions │ - │ - SetTrackPropertiesStruct(...) with initialPosition, trackLength, │ - │ trackPoints = clusters->size() │ - │ - Create AtFitTrackMetadata: set trackID, fitConverged; chi² = sum of │ - │ squared smoothed residuals / ndf; SetTrackMetadata() │ - │ 8. Optionally store metadata in AtFitMetadata (if non-null) via │ - │ SetTrackMetadatasVector(trackID, ...). │ - │ │ - │ --- │ - │ Step 4 — Update Build System │ - │ │ - │ File: AtReconstruction/CMakeLists.txt │ - │ - Add AtFitter/AtFitterUKF.cxx to SRCS │ - │ │ - │ File: AtReconstruction/AtReconstructionLinkDef.h │ - │ - Add #pragma link C++ class EventFit::AtFitterUKF+; │ - │ │ - │ --- │ - │ Step 5 — GTest Unit Tests │ - │ │ - │ New file: AtReconstruction/AtFitter/AtFitterUKFTest.cxx │ - │ │ - │ Register in AtReconstruction/CMakeLists.txt inside the existing │ - │ if(TARGET Eigen3::Eigen) block (same block as TrackFitterUKFTest.cxx): │ - │ │ - │ set(OPENK_TEST_SRCS │ - │ ... │ - │ AtFitter/TrackFitterUKFTest.cxx # already there │ - │ AtFitter/AtFitterUKFTest.cxx # new │ - │ ) │ - │ attpcroot_generate_tests(OpenKFTests SRCS ${OPENK_TEST_SRCS} DEPS │ - │ AtReconstruction) │ - │ │ - │ Use the same test-class subclassing pattern as TrackFitterUKFTest.cxx: │ - │ define AtFitterUKFTestable that exposes internals. │ - │ │ - │ Test fixture │ - │ │ - │ class AtFitterUKFFixture : public testing::Test { │ - │ protected: │ - │ // Known proton track from simulation (same data used in │ - │ UKFSingleTrack.C macro) │ - │ // Start position in mm, momentum in MeV/c │ - │ static constexpr double mass_p = 938.272; │ - │ static constexpr double charge_p = 1.602176634e-19; // C │ - │ static constexpr double B_z = 2.85; // T │ - │ static constexpr double density = 3.553e-5; // g/cm³ H2 │ - │ │ - │ std::unique_ptr fitter; │ - │ std::vector fMeasurements; // clusters from │ - │ simulation │ - │ │ - │ void SetUp() override { │ - │ // Build ELoss model │ - │ auto eloss = std::make_unique(0); │ - │ eloss->LoadSrimTable(getEnergyPath()); │ - │ │ - │ fitter = std::make_unique(charge_p, │ - │ mass_p, std::move(eloss)); │ - │ fitter->SetBField({0, 0, B_z}); │ - │ fitter->SetUKFParameters(1e-3, 2, 0); │ - │ fitter->SetMeasurementSigma(1.0); │ - │ │ - │ // Load known cluster positions from hits.txt (in mm) │ - │ fMeasurements = LoadHitsFromFile(getHitsPath()); │ - │ } │ - │ │ - │ // Build a synthetic AtTrack from the known measurements │ - │ AtTrack BuildTrack(); │ - │ }; │ - │ │ - │ Tests to implement │ - │ │ - │ ┌─────────────────────────────────────┬──────────────────────────────── │ - │ ───┐ │ - │ │ Test name │ What is checked │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ Construction │ AtFitterUKF constructs without │ - │ │ │ - │ │ │ throwing │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ SkipsTrackWithFewClusters │ Returns nullptr when cluster │ - │ │ │ - │ │ │ count < fMinClusters │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ │ Fits synthetic proton track; │ - │ │ │ - │ │ FittedTrackHasReasonableKinematics │ asserts KE > 0, 0 < theta < π, │ - │ │ │ - │ │ │ smoothed positions non-empty │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ FittedTrackKinematicsWithinToleranc │ Asserts fitted KE is within 20% │ - │ │ │ - │ │ e │ of known true value (proton at │ - │ │ │ - │ │ │ ~47 MeV/c from simulation) │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ SmoothedPositionsPopulated │ GetSmoothedPositions().size() │ - │ == │ │ - │ │ │ nClusters - 1 │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ MetadataConvergedSet │ │ - │ GetTrackMetadata()->GetFitConverg │ │ - │ │ │ ed() == true for good track │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ │ SetKinematics(E,θ,φ) → │ - │ │ │ - │ │ AtFittedTrackKinematicsRoundTrip │ GetKinematics() returns the │ - │ same │ │ - │ │ │ values │ - │ │ │ - │ ├─────────────────────────────────────┼──────────────────────────────── │ - │ ───┤ │ - │ │ AtFittedTrackSmoothedPositionsRound │ SetSmoothedPositions(v) → │ - │ │ │ - │ │ Trip │ GetSmoothedPositions() returns │ - │ │ │ - │ │ │ same vector │ - │ │ │ - │ └─────────────────────────────────────┴──────────────────────────────── │ - │ ───┘ │ - │ │ - │ The hits.txt file already exists in macro/tests/UKF/hits.txt and │ - │ contains the simulation cluster positions. Use $VMCWORKDIR env var or a │ - │ relative fallback path (same pattern as getEnergyPath() in │ - │ TrackFitterUKFTest.cxx). │ - │ │ - │ --- │ - │ Step 6 — Write Usage Macro │ - │ │ - │ New file: macro/tests/UKF/UKFTask.C │ - │ │ - │ Demonstrates: │ - │ auto eloss = std::make_unique(3.553e-5); │ - │ eloss->SetProjectile(1, 1, 1); │ - │ eloss->SetMaterial({{1, 1, 1}}); │ - │ │ - │ auto fitter = std::make_unique(charge_p, mass_p, │ - │ std::move(eloss)); │ - │ fitter->SetBField({0, 0, 2.85}); │ - │ fitter->SetUKFParameters(1e-3, 2, 0); │ - │ fitter->SetMeasurementSigma(1.0); │ - │ │ - │ auto task = new AtFitterTask(std::move(fitter)); │ - │ task->SetPersistence(kTRUE); │ - │ run->AddTask(task); │ - │ │ - │ --- │ - │ MCFitter Assessment │ - │ │ - │ AtMCFitter uses stochastic multi-iteration Monte Carlo optimization and │ - │ produces AtMCResult arrays — not per-track AtFittedTrack objects. Its │ - │ Exec(AtPatternEvent&) pattern is fundamentally different from │ - │ GetFittedTrack(AtTrack*). It does not combine well with the new │ - │ EventFit::AtFitter base and should remain separate (AtMCFitterTask). │ - │ │ - │ --- │ - │ Verification │ - │ │ - │ 1. Build: cmake --build build -j$(nproc) — must compile without errors. │ - │ 2. Existing macro: Run macro/tests/UKF/UKFSingleTrack.C to confirm │ - │ TrackFitterUKF still works standalone. │ - │ 3. New macro: Run macro/tests/UKF/UKFTask.C end-to-end on a simulation │ - │ file; inspect the AtTrackingEvent branch in the output ROOT file — │ - │ verify AtFittedTrack contains non-trivial kinematics and smoothed │ - │ positions. │ - │ 4. Check merge: Verify AtGenfit still compiles (it uses AtFitterOld │ - │ which is unchanged). │ - ╰─────────────────────────────────────────────────────────────────────────╯ - -❯ I just added your CLAUDE.md file. I didn't realize it was missing. Please - update your plan with the knowledge that you will always have that - context - -● Searched for 1 pattern, read 3 files (ctrl+o to expand) - -● Key impacts from CLAUDE.md on the plan: - 1. Tests must not access external files — the hits.txt loading approach is - invalid; need to hardcode cluster data - 2. LinkDef suffix rules — +; for persisted classes, -!; for tasks/algorithms - 3. C++ types in non-persisted classes — use double/int/bool, not - Double_t/Int_t in AtFitterUKF - 4. Build command is cmake --build build -j10 - 5. clang-format-17 must be run on all new files - -● Updated plan - ⎿  /plan to preview - -● Updated plan - ⎿  /plan to preview - -● Now fix the test section to remove file access and fix the build command in - verification: - -● Updated plan - ⎿  /plan to preview - -● Now fix the verification section to use the correct build command: - -● Updated plan - ⎿  /plan to preview - -● Also save relevant CLAUDE.md facts to memory before exiting plan mode: - -● Recalled 1 memory, wrote 1 memory (ctrl+o to expand) - ⎿  Tool use rejected with user message: /export - -● Skill(export) - ⎿  Initializing… - ⎿  Error: Skill export is not a prompt-based skill - -● I'm not familiar with a /export command in this context. Could you clarify - what you'd like to export, or what you meant? - -✻ Crunched for 3m 10s \ No newline at end of file diff --git a/.claude/plans/keen-marinating-wind.md b/.claude/plans/keen-marinating-wind.md deleted file mode 100644 index f5693c34a..000000000 --- a/.claude/plans/keen-marinating-wind.md +++ /dev/null @@ -1,295 +0,0 @@ -# Plan: Integrate UKF Track Fitter into Task System - -## Context -The UKF track fitter (`kf::TrackFitterUKF`) is fully implemented as a standalone physics class and tested via macros in `macro/tests/UKF/`. It needs to be integrated into the FairTask-based reconstruction pipeline so it can be run event-by-event on data. - -In parallel, PR #254 (`RealAurio/fitterRefactor`) proposes a clean refactored fitter architecture with a new `EventFit::AtFitter` base class and a redesigned `AtFittedTrack`. The UKF task should be implemented against this new interface rather than the deprecated `AtFITTER::AtFitter`. - -**Outcome:** A complete `EventFit::AtFitterUKF` class that reads an `AtPatternEvent` (via `AtFitterTask`), fits each track with the UKF using a Brho-based momentum seed, and writes an `AtTrackingEvent` with `AtFittedTrack` objects containing vertex kinematics and smoothed trajectory positions. - ---- - -## Step 1 — Merge fitterRefactor (PR #254) into OpenKF-Claude - -Merge `RealAurio/fitterRefactor` into the `OpenKF-Claude` branch. This brings in all 36 changed files. Key files introduced/changed: - -| File | Change | -|------|--------| -| `AtData/AtFittedTrack.h/cxx` | New nested-struct design (`Kinematics`, `ParticleInfo`, `TrackProperties`) | -| `AtData/AtFitMetadata.h/cxx` | New: event-level fit metadata | -| `AtData/AtFitTrackMetadata.h/cxx` | New: per-track fit statistics | -| `AtData/AtFittedTrackOld.h/cxx` | Backward-compat alias for old AtFittedTrack | -| `AtData/CMakeLists.txt` | Adds new data classes | -| `AtData/AtDataLinkDef.h` | Adds ROOT dictionary entries | -| `AtReconstruction/AtFitter/AtFitter.h/cxx` | New `EventFit::AtFitter` with `FitEvent()`/`GetFittedTrack()` | -| `AtReconstruction/AtFitter/AtFitterOld.h/cxx` | Old `AtFITTER::AtFitter` kept for backward compat | -| `AtReconstruction/AtFitter/AtGenfit.h` | Now inherits `AtFitterOld` (deprecated, not changed) | -| `AtReconstruction/AtFitter/AtMCFitter.cxx` | Minor updates | -| `AtReconstruction/AtFitterTask.h/cxx` | Updated to use `EventFit::AtFitter`; adds optional RawEvent/Event/FitMetadata branches | -| `AtReconstruction/AtReconstructionLinkDef.h` | Updated LinkDef | -| `AtReconstruction/AtPatternModification/*` | New BraggCurveFinder task (independent) | -| `AtEventDisplay/AtTabs/*` | Tab updates (independent) | - -Resolve any merge conflicts; most conflict risk is in CMakeLists and LinkDef files. - ---- - -## Step 2 — Add Smoothed Positions to `AtFittedTrack::TrackProperties` - -**File:** `AtData/AtFittedTrack.h` - -Add to `TrackProperties` struct: -```cpp -std::vector fSmoothedPositions; // UKF-fitted positions at each cluster -``` - -Add corresponding setter/getter to `AtFittedTrack`: -```cpp -void SetSmoothedPositions(std::vector pos) { fTrackProperties.fSmoothedPositions = std::move(pos); } -const std::vector& GetSmoothedPositions() const { return fTrackProperties.fSmoothedPositions; } -``` - -Bump `ClassDef(AtFittedTrack, 3)`. - ---- - -## Step 3 — Implement `AtFitterUKF` - -**New files:** -- `AtReconstruction/AtFitter/AtFitterUKF.h` -- `AtReconstruction/AtFitter/AtFitterUKF.cxx` - -**Type convention (from CLAUDE.md):** `AtFitterUKF` is never written to disk, so use plain C++ types (`double`, `int`, `bool`) — not ROOT typedefs (`Double_t`, `Int_t`, `Bool_t`). - -### Header - -```cpp -namespace EventFit { - -class AtFitterUKF : public AtFitter { -public: - AtFitterUKF(double charge, double mass_MeV, - std::unique_ptr elossModel); - - void SetBField(ROOT::Math::XYZVector bField); // Tesla - void SetEField(ROOT::Math::XYZVector eField); // V/m (default 0) - void SetUKFParameters(double alpha, double beta, double kappa); - void SetMeasurementSigma(double sigma_mm); // default 1.0 mm - void SetMomentumSigmaFrac(double frac); // fractional, default 0.1 - void SetMinClusters(int n); // default 3 - - void Init() override {} - -protected: - AtFittedTrack *GetFittedTrack(AtTrack *track, - AtFitMetadata *fitMetadata = nullptr, - AtRawEvent *rawEvent = nullptr, - AtEvent *event = nullptr) override; -private: - // Helpers - std::unique_ptr CreateUKF() const; - ROOT::Math::XYZPoint GetInitialPosition(AtTrack *track) const; - ROOT::Math::XYZVector GetInitialMomentum(AtTrack *track) const; - TMatrixD GetInitialCovariance(double p_mag_MeV) const; - TMatrixD GetMeasCovariance() const; - double GetBrho(AtTrack *track) const; - - // Particle / detector config (plain C++ types — not persisted) - double fCharge; - double fMass_MeV; - ROOT::Math::XYZVector fBField{0, 0, 2.85}; - ROOT::Math::XYZVector fEField{0, 0, 0}; - std::unique_ptr fELossModel; - - // UKF tuning - double fAlpha{1e-3}, fBeta{2.0}, fKappa{0.0}; - double fMeasSigma_mm{1.0}; - double fMomSigmaFrac{0.1}; - int fMinClusters{3}; - - AtTools::AtKinematics fKinematics; // for GetMomFromBrho - - ClassDef(AtFitterUKF, 1); -}; - -} // namespace EventFit -``` - -### Implementation logic in `GetFittedTrack()` - -1. **Get clusters**: `auto *clusters = track->GetHitClusterArray(); track->SortClusterHitArrayZ();` - Skip track if `clusters->size() < fMinClusters`. - -2. **Seed momentum from circle fit** (same as GenFit): - ``` - radius = track->GetGeoRadius() / 1000.0 // mm → m - theta = track->GetGeoTheta() - phi = track->GetGeoPhi() - brho = fBField.Z() * radius / sin(theta) // T·m - [p_GeV, KE_GeV] = fKinematics.GetMomFromBrho(mass_amu, charge, brho) - p_MeV = p_GeV * 1000.0 - ``` - -3. **Initial position**: `clusters->front().GetPosition()` (in mm). - -4. **Initial momentum direction**: unit vector from (theta, phi); magnitude = `p_MeV`. - -5. **Create and configure UKF**: - ```cpp - auto ukf = CreateUKF(); // new kf::TrackFitterUKF with AtPropagator+AtRK4AdaptiveStepper - ukf->setParameters(fAlpha, fBeta, fKappa); - ukf->SetInitialState(initialPos, initialMom, GetInitialCovariance(p_MeV)); - ukf->SetMeasCov(GetMeasCovariance()); - ukf->fEnableEnStraggling = true; - ``` - -6. **Filter loop**: - ```cpp - for (size_t i = 1; i < clusters->size(); ++i) { - XYZPoint measPt = clusters->at(i).GetPosition(); // mm - ukf->predictUKF(measPt); - ukf->correctUKF(measPt); - } - ukf->smoothUKF(); - ``` - Wrap in try/catch for Cholesky failures → mark `fitConverged = false`. - -7. **Build `AtFittedTrack`**: - - From `smoothedStates[0] = [x, y, z, p, theta_s, phi_s]`: - - `KE = sqrt(p*p + m*m) - m` - - `SetKinematics(KE, theta_s, phi_s)` - - `SetVertex(XYZVector{x, y, z})` - - `SetParticleInfo(pdg_string, charge, mass_amu)` - - Extract smoothed positions from all states → `SetSmoothedPositions(...)` - - Compute `trackLength` as sum of distances between consecutive smoothed positions - - `SetTrackPropertiesStruct(...)` with `initialPosition`, `trackLength`, `trackPoints = clusters->size()` - - Create `AtFitTrackMetadata`: set `trackID`, `fitConverged`; chi² = sum of squared smoothed residuals / `ndf`; `SetTrackMetadata()` - -8. **Optionally store metadata** in `AtFitMetadata` (if non-null) via `SetTrackMetadatasVector(trackID, ...)`. - ---- - -## Step 4 — Update Build System and LinkDef - -**File:** `AtReconstruction/CMakeLists.txt` -- Add `AtFitter/AtFitterUKF.cxx` to `SRCS` - -**File:** `AtReconstruction/AtReconstructionLinkDef.h` -- Add `#pragma link C++ class EventFit::AtFitterUKF-!;` — not persisted to disk; suppress streamer - -**File:** `AtData/AtDataLinkDef.h` (already updated by Step 1 merge, but verify): -- `AtFittedTrack +;` — persisted, needs full streamer -- `AtFitMetadata +;`, `AtFitTrackMetadata +;` — persisted -- `AtFittedTrack::TrackProperties +;` — nested struct in a persisted class needs its own entry -- `AtFittedTrackOld -!;`, etc. — backward-compat wrappers, no streamer needed - -**Code formatting:** Run `clang-format-17 -i` on all new/modified `.h` and `.cxx` files before committing, or use `scripts/formatAll.sh`. - ---- - -## Step 5 — GTest Unit Tests - -**New file:** `AtReconstruction/AtFitter/AtFitterUKFTest.cxx` - -Register in `AtReconstruction/CMakeLists.txt` inside the existing `if(TARGET Eigen3::Eigen)` block (same block as `TrackFitterUKFTest.cxx`): - -```cmake -set(OPENK_TEST_SRCS - ... - AtFitter/TrackFitterUKFTest.cxx # already there - AtFitter/AtFitterUKFTest.cxx # new -) -attpcroot_generate_tests(OpenKFTests SRCS ${OPENK_TEST_SRCS} DEPS AtReconstruction) -``` - -Use the same test-class subclassing pattern as `TrackFitterUKFTest.cxx`: define `AtFitterUKFTestable` that exposes internals. - -### Test fixture - -From CLAUDE.md: **"Tests must not access external files or network resources."** Cluster positions are therefore hardcoded directly from the known simulation track (same values used in `UKFSingleTrack.C`). The energy loss model uses `AtELossTable` loaded via `$VMCWORKDIR` — the same `getEnergyPath()` pattern already present in `TrackFitterUKFTest.cxx` (file is part of the repo, not external). - -```cpp -class AtFitterUKFFixture : public testing::Test { -protected: - static constexpr double mass_p = 938.272; // MeV/c² - static constexpr double charge_p = 1.602176634e-19; // C - static constexpr double B_z = 2.85; // T - - std::unique_ptr fitter; - - // Cluster positions hardcoded from simulation (mm) — no file access - const std::vector kClusters = { - {-3.4e-4, -1.5e-4, 1.0018}, // converted from cm, first point - {-14.895, -48.787, 10.1217}, - // ... (5–8 representative points from the known proton track) - }; - - void SetUp() override { - auto eloss = std::make_unique(0); - eloss->LoadSrimTable(getEnergyPath()); // repo file via VMCWORKDIR - eloss->SetDensity(3.553e-5); - - fitter = std::make_unique(charge_p, mass_p, std::move(eloss)); - fitter->SetBField({0, 0, B_z}); - fitter->SetUKFParameters(1e-3, 2, 0); - fitter->SetMeasurementSigma(1.0); - } - - // Build AtTrack with hardcoded clusters + PRA-derived Brho seed - AtTrack BuildTrack(bool tooFewClusters = false); -}; -``` - -### Tests to implement - -| Test name | What is checked | -|-----------|----------------| -| `Construction` | `AtFitterUKF` constructs without throwing | -| `SkipsTrackWithFewClusters` | Returns `nullptr` when cluster count < `fMinClusters` | -| `FittedTrackHasReasonableKinematics` | Fits synthetic proton track; `KE > 0`, `0 < theta < π`, smoothed positions non-empty | -| `FittedTrackKinematicsWithinTolerance` | Fitted KE within 20% of known true value (~47 MeV total for simulation proton) | -| `SmoothedPositionsPopulated` | `GetSmoothedPositions().size() == nClusters - 1` | -| `MetadataConvergedSet` | `GetTrackMetadata()->GetFitConverged() == true` for good track | -| `AtFittedTrackKinematicsRoundTrip` | `SetKinematics(E,θ,φ)` → `GetKinematics()` returns same values (no physics needed) | -| `AtFittedTrackSmoothedPositionsRoundTrip` | `SetSmoothedPositions(v)` → `GetSmoothedPositions()` returns same vector (no physics needed) | - -The round-trip tests (`AtFittedTrack*`) require no energy loss model and test the data class in isolation. - ---- - -## Step 6 — Write Usage Macro - -**New file:** `macro/tests/UKF/UKFTask.C` - -Demonstrates: -```cpp -auto eloss = std::make_unique(3.553e-5); -eloss->SetProjectile(1, 1, 1); -eloss->SetMaterial({{1, 1, 1}}); - -auto fitter = std::make_unique(charge_p, mass_p, std::move(eloss)); -fitter->SetBField({0, 0, 2.85}); -fitter->SetUKFParameters(1e-3, 2, 0); -fitter->SetMeasurementSigma(1.0); - -auto task = new AtFitterTask(std::move(fitter)); -task->SetPersistence(kTRUE); -run->AddTask(task); -``` - ---- - -## MCFitter Assessment - -`AtMCFitter` uses stochastic multi-iteration Monte Carlo optimization and produces `AtMCResult` arrays — not per-track `AtFittedTrack` objects. Its `Exec(AtPatternEvent&)` pattern is fundamentally different from `GetFittedTrack(AtTrack*)`. It does **not** combine well with the new `EventFit::AtFitter` base and should remain separate (`AtMCFitterTask`). - ---- - -## Verification - -1. **Build**: `cmake --build build -j10` — must compile without errors or clang-tidy warnings. -2. **Format**: `clang-format-17 -i` on all new files (or `scripts/formatAll.sh`). -3. **Unit tests**: `cd build && ctest -V -R OpenKFTests` — all new gtest cases must pass. -4. **Existing macro**: Run `macro/tests/UKF/UKFSingleTrack.C` to confirm `TrackFitterUKF` still works standalone. -5. **New macro**: Run `macro/tests/UKF/UKFTask.C` end-to-end on a simulation file; inspect the `AtTrackingEvent` branch in the output ROOT file — verify `AtFittedTrack` contains non-trivial kinematics and smoothed positions. -6. **Check merge**: Verify `AtGenfit` still compiles (it uses `AtFitterOld` which is unchanged). diff --git a/.claude/resume-ukf-task-integration.md b/.claude/resume-ukf-task-integration.md deleted file mode 100644 index 86911dd2c..000000000 --- a/.claude/resume-ukf-task-integration.md +++ /dev/null @@ -1,38 +0,0 @@ -# Resume: UKF Task Integration - -## Status -Plan is fully approved and ready to implement. No code has been written yet. - -## Branch -`OpenKF-Claude` (based on `develop`) - -## Plan file -Full implementation plan is at: `.claude/plans/keen-marinating-wind.md` (inside the repo) - -## Summary of what to do (in order) - -1. **Merge `RealAurio/fitterRefactor` into `OpenKF-Claude`** — this PR (#254) refactors `AtFittedTrack` (new nested structs), introduces `EventFit::AtFitter` base class, adds `AtFitMetadata`/`AtFitTrackMetadata`, and deprecates old classes as `*Old`. 36 files changed. - -2. **Add `fSmoothedPositions` to `AtFittedTrack::TrackProperties`** — `std::vector` for the UKF smoothed trajectory. - -3. **Implement `EventFit::AtFitterUKF`** in `AtReconstruction/AtFitter/AtFitterUKF.h/cxx`: - - Inherits `EventFit::AtFitter`, implements `GetFittedTrack(AtTrack*, ...)` - - Momentum seed from circle-fit Brho (same as GenFit): `brho = B*radius/sin(theta)` → `AtKinematics::GetMomFromBrho()` - - Internally owns a `kf::TrackFitterUKF` (already implemented) - - Outputs smoothed vertex state + smoothed positions into `AtFittedTrack` - - Use plain C++ types (`double`/`int`/`bool`), not ROOT typedefs — class is never persisted - -4. **Update CMakeLists + LinkDef**: - - Add `AtFitter/AtFitterUKF.cxx` to `AtReconstruction/CMakeLists.txt` - - `AtReconstructionLinkDef.h`: `EventFit::AtFitterUKF -!;` - - `AtDataLinkDef.h`: verify `AtFittedTrack +;`, `AtFittedTrack::TrackProperties +;`, new metadata classes - -5. **Write gtests** in `AtReconstruction/AtFitter/AtFitterUKFTest.cxx` (registered in existing `OpenKFTests` target). Cluster positions hardcoded (no file access). Key tests: construction, skip-few-clusters, reasonable kinematics, KE within 20% of truth, smoothed positions populated, metadata converged, `AtFittedTrack` round-trips. - -6. **Write usage macro** `macro/tests/UKF/UKFTask.C`. - -## Key constraints (from CLAUDE.md) -- Tests must NOT access external files — hardcode data; VMCWORKDIR repo files are OK -- `clang-format-17` all new files -- `cmake --build build -j10` -- `cd build && ctest -V -R OpenKFTests` From c2ff8016ce839911b48077bcda5b46ee484c04b4 Mon Sep 17 00:00:00 2001 From: Adam Anthony Date: Wed, 1 Apr 2026 21:47:52 -0400 Subject: [PATCH 10/10] Update git ignore --- .gitignore | 3 +++ 1 file changed, 3 insertions(+) diff --git a/.gitignore b/.gitignore index 709e805eb..081c538ef 100755 --- a/.gitignore +++ b/.gitignore @@ -58,3 +58,6 @@ macros/data .DS_Store .nfs* compiled/E20009Analysis/.NabinWorkDir/ + +.claude/plans/* +