mx/Jungfraujoch
3
2 Fork 0
Code Issues Pull Requests Actions Packages Projects Releases 89 Wiki Activity

rugnux: panel study log of held-out residual per panel and shell (RUGNUX_PANEL_STUDY)

Browse Source 
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01UqXYdJKPUFLoeAchjzv3x4
This commit is contained in:
leonarski_fandClaude Opus 5.5 committed 2026-11-03 07:15:21 +01:00
1 parent a82037a41b
commit ffc990a330
1 file changed
+79 -11
Download Patch File Download Diff File Expand all files Collapse all files
rugnux/RugnuxStills.cpp
+79 -11
View File
@@ -116,6 +116,59 @@ GeomRefineFrame ToRefineFrame(const DataMessage &msg, int ordinal) {
return f;
}
// One held-out fit: the geometry fitted on one half, and its residuals on the other half's spots.
struct HeldOutFit {
double rms = std::numeric_limits<double>::infinity();
std::vector<ModuleShift> modules; // the modules as fitted on the fitting half
GeometryRefinerResult eval; // the frozen geometry on the other half
DiffractionGeometry geom;
};
// PANEL STUDY (RUGNUX_PANEL_STUDY): the held-out residual of each fold per module and per resolution shell,
// and each module's shift on the two halves, to the log under "PanelStudy".
void LogPanelStudy(Logger &logger, const std::string &label, const std::vector<GeomRefineFrame> half[2],
const HeldOutFit fits[2]) {
constexpr double CAP_PX = 2.0;
const std::vector<double> edges = {1e9, 4.0, 3.0, 2.5, 2.2, 2.0, 1.8, 1.6, 1.4, 0.0};
const size_t nm = fits[0].modules.size();
std::vector<double> m_sq(nm, 0), s_sq(edges.size(), 0);
std::vector<int> m_n(nm, 0), m_hit(nm, 0), s_n(edges.size(), 0), s_hit(edges.size(), 0);
for (int fit = 0; fit < 2; fit++) {
const auto &e = fits[fit].eval;
const auto &frames = half[1 - fit];
if (e.predicted_px.size() != frames.size())
continue;
for (size_t i = 0; i < frames.size(); i++)
for (size_t j = 0; j < frames[i].spots.size() && j < e.predicted_px[i].size(); j++) {
const auto &sp = frames[i].spots[j];
const auto &pr = e.predicted_px[i][j];
double r = CAP_PX;
if (std::isfinite(pr[0]) && std::isfinite(pr[1]))
r = std::min(CAP_PX, static_cast<double>(std::hypot(pr[0] - sp.x, pr[1] - sp.y)));
const double d = 1.0 / std::max(1e-6f, fits[fit].geom.DetectorToRecip(sp.x, sp.y).Length());
size_t sh = 1;
while (sh + 1 < edges.size() && d < edges[sh])
sh++;
s_sq[sh] += r * r; s_n[sh]++; s_hit[sh] += r < 1.0;
if (sp.module >= 0 && static_cast<size_t>(sp.module) < nm) {
m_sq[sp.module] += r * r; m_n[sp.module]++; m_hit[sp.module] += r < 1.0;
}
}
}
for (size_t sh = 1; sh < edges.size(); sh++)
if (s_n[sh] > 0)
logger.Info("PanelStudy {} shell {:.2f} {:.2f} n {} rms {:.4f} hit1 {:.4f}", label, edges[sh - 1], edges[sh],
s_n[sh], std::sqrt(s_sq[sh] / s_n[sh]), s_hit[sh] / static_cast<double>(s_n[sh]));
for (size_t m = 0; m < nm; m++)
if (m_n[m] > 0)
logger.Info("PanelStudy {} module {} x0 {:.0f} y0 {:.0f} n {} rms {:.4f} hit1 {:.4f} "
"fold0 {:+.3f} {:+.3f} {} fold1 {:+.3f} {:+.3f} {}", label, m,
fits[0].modules[m].x0, fits[0].modules[m].y0, m_n[m], std::sqrt(m_sq[m] / m_n[m]),
m_hit[m] / static_cast<double>(m_n[m]),
fits[0].modules[m].dx, fits[0].modules[m].dy, fits[0].modules[m].spots,
fits[1].modules[m].dx, fits[1].modules[m].dy, fits[1].modules[m].spots);
}
bool SameGeometry(const DiffractionGeometry &a, const DiffractionGeometry &b) {
return a.GetBeamX_pxl() == b.GetBeamX_pxl() && a.GetBeamY_pxl() == b.GetBeamY_pxl()
&& a.GetDetectorDistance_mm() == b.GetDetectorDistance_mm() && a.GetPixelSize_mm() == b.GetPixelSize_mm()
@@ -506,32 +559,47 @@ void Rugnux::RefineStillsGeometry(int start_image, int end_image, int images_to_
});
// Held-out residual of the geometry fitted on half `fit`, scored on every spot of the other half.
auto held_out = [&](const std::vector<ModuleShift> &modules, bool refine, int fit) {
auto held_out = [&](const std::vector<ModuleShift> &modules, bool refine, int fit, int min_spots) {
HeldOutFit h;
GeometryRefinerSettings sf = gr_settings;
sf.num_threads = std::max(1, config_.nthreads / 8);
sf.modules = modules;
sf.refine_modules = refine;
sf.min_module_spots = min_spots;
const auto f = RefineGlobalGeometry(nominal, *cell, half[fit], sf);
if (!f.ok)
return std::numeric_limits<double>::infinity();
return h;
DiffractionGeometry g = nominal;
g.BeamX_pxl(f.beam_x_px).BeamY_pxl(f.beam_y_px).DetectorDistance_mm(f.distance_mm);
GeometryRefinerSettings ev = gr_settings;
ev.num_threads = sf.num_threads;
ev.refine_beam_distance = false;
ev.modules = f.modules;
const auto e = RefineGlobalGeometry(g, f.cell, half[1 - fit], ev);
if (!e.ok)
return std::numeric_limits<double>::infinity();
return e.capped_rms_px;
h.modules = f.modules;
h.geom = g;
h.eval = RefineGlobalGeometry(g, f.cell, half[1 - fit], ev);
if (h.eval.ok)
h.rms = h.eval.capped_rms_px;
return h;
};
std::future<double> held_fits[2], moved_fits[2];
const int min_spots = gr_settings.min_module_spots;
std::future<HeldOutFit> held_fits[2], moved_fits[2];
for (int fit = 0; fit < 2; fit++) {
held_fits[fit] = std::async(std::launch::async, held_out, std::cref(centres), false, fit);
moved_fits[fit] = std::async(std::launch::async, held_out, std::cref(start), true, fit);
held_fits[fit] = std::async(std::launch::async, held_out, std::cref(centres), false, fit, min_spots);
moved_fits[fit] = std::async(std::launch::async, held_out, std::cref(start), true, fit, min_spots);
}
HeldOutFit held_f[2] = {held_fits[0].get(), held_fits[1].get()};
HeldOutFit moved_f[2] = {moved_fits[0].get(), moved_fits[1].get()};
const double held = (held_f[0].rms + held_f[1].rms) / 2.0;
const double moved = (moved_f[0].rms + moved_f[1].rms) / 2.0;
if (std::getenv("RUGNUX_PANEL_STUDY")) {
LogPanelStudy(logger, "given", half, held_f);
LogPanelStudy(logger, "placed", half, moved_f);
HeldOutFit low_f[2] = {held_out(start, true, 0, 100), held_out(start, true, 1, 100)};
LogPanelStudy(logger, "placed100", half, low_f);
logger.Info("PanelStudy summary given {:.4f} placed {:.4f} placed100 {:.4f}", held, moved,
(low_f[0].rms + low_f[1].rms) / 2.0);
}
const double held = (held_fits[0].get() + held_fits[1].get()) / 2.0;
const double moved = (moved_fits[0].get() + moved_fits[1].get()) / 2.0;
const GeometryRefinerResult rm = placement.get();
placement_mark.Stop();
logger.Info("Module placement: {} modules; held-out spot residual (capped rms) {:.3f} px as given, "