diff --git a/common/Reflection.h b/common/Reflection.h index 1b8b5caec..1ab59fc19 100644 --- a/common/Reflection.h +++ b/common/Reflection.h @@ -62,6 +62,9 @@ struct Reflection { // intensity at all - its brightest pixels are missing - and a rotation merge drops the whole // rocking event it belongs to, as XDS does with an overloaded reflection. bool overloaded = false; + // A still's overloaded reflection kept as a lower bound instead (BraggIntegrationEngine::Finalize); + // the run may still hand it back to the overload drop (RugnuxScaleMerge). + bool overload_kept = false; }; // One full reflection of a rotation sweep as the merge used it - its partials summed into one diff --git a/image_analysis/bragg_integration/BraggIntegrationEngine.cpp b/image_analysis/bragg_integration/BraggIntegrationEngine.cpp index b9cc10859..8ed7b99a6 100644 --- a/image_analysis/bragg_integration/BraggIntegrationEngine.cpp +++ b/image_analysis/bragg_integration/BraggIntegrationEngine.cpp @@ -288,6 +288,7 @@ std::vector BraggIntegrationEngine::Finalize(const std::vector(fr.saturated * (saturation_limit - fr.bkg)); refl.sigma = std::sqrt(fr.var_bkg + std::max(refl.I, 1.0f)); refl.overloaded = false; + refl.overload_kept = true; ++counts.overloads_kept; } refl.observed = true; diff --git a/rugnux/RugnuxScaleMerge.cpp b/rugnux/RugnuxScaleMerge.cpp index 815b0edec..b3e876321 100644 --- a/rugnux/RugnuxScaleMerge.cpp +++ b/rugnux/RugnuxScaleMerge.cpp @@ -549,6 +549,45 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) { size_t ice_flagged = 0, ice_total = 0; { auto &outcomes = indexer->GetIntegrationOutcome(); + // Saturation clips the strongest spots, and those are at low resolution: kept overloads of a + // still are a measurement only where they concentrate there. Where they do not, the saturated + // pixels are something else (hot or flagged pixels across the pattern), and an observation + // with such a pixel set to the saturation limit is not a lower bound but an outlier - so they + // go back to the overload drop. A pragmatic gate until saturation is told from bad pixels per + // spot. + { + std::atomic k_all{0}, k_low{0}, n_all{0}, n_low{0}; + ParallelFor(static_cast(outcomes.size()), config_.nthreads, [&](int o) { + int64_t all = 0, low = 0, kept = 0, kept_low = 0; + for (const auto &r : outcomes[o].reflections) { + const bool low_res = r.d >= OVERLOAD_LOW_RES_A; + ++all; + low += low_res; + if (r.overload_kept) { + ++kept; + kept_low += low_res; + } + } + n_all += all; n_low += low; k_all += kept; k_low += kept_low; + }); + if (k_all > 0 && n_low > 0) { + const double frac_all = static_cast(k_all) / static_cast(n_all); + const double frac_low = static_cast(k_low) / static_cast(n_low); + if (frac_low < 1.5 * frac_all) { + for (auto &oc : outcomes) + for (auto &r : oc.reflections) + if (r.overload_kept) { + r.overloaded = true; + r.overload_kept = false; + } + result.overloads_kept = 0; + logger.Info("Overloads: the {} kept as a lower bound are not concentrated at low resolution " + "({:.3f}% of the reflections coarser than {:.0f} A, {:.3f}% of all) - not " + "saturated spots, so they are dropped as overloads", + k_all.load(), 100.0 * frac_low, OVERLOAD_LOW_RES_A, 100.0 * frac_all); + } + } + } std::atomic n_all{0}, n_low{0}, n_over{0}, n_over_low{0}, n_ice{0}; ParallelFor(static_cast(outcomes.size()), config_.nthreads, [&](int o) { int64_t all = 0, low = 0, over = 0, over_low = 0, on_ice = 0;