diff --git a/image_analysis/scale_merge/SearchSpaceGroup.cpp b/image_analysis/scale_merge/SearchSpaceGroup.cpp index bb3f9c1c..935df8ec 100644 --- a/image_analysis/scale_merge/SearchSpaceGroup.cpp +++ b/image_analysis/scale_merge/SearchSpaceGroup.cpp @@ -219,14 +219,41 @@ SearchSpaceGroupResult SearchSpaceGroup( const bool in_range = finite && (opt.d_min_limit_A <= 0 || r.d >= opt.d_min_limit_A); IoverSigma[i] = finite ? r.I / r.sigma : 0.0; pass_absence[i] = in_range; - // The correlation stage uses only genuinely-present reflections. Near-zero (systematically - // absent) reflections would otherwise form a second cluster at the origin and fake a high - // correlation for false operators - fatal on centered lattices, where half the reflections - // are extinct. - pass_cc[i] = in_range && IoverSigma[i] >= opt.present_i_over_sigma && - (opt.min_i_over_sigma <= 0 || IoverSigma[i] >= opt.min_i_over_sigma); } + // The merge floors the merged sigma at b|I| (Merge.h SigmaWithSystematicFloor) so that ISa = 1/b + // stays the asymptotic I/sigma ceiling: no reflection in the merge can read above it. A fixed + // I/sigma cut is therefore not a per-reflection test but a switch on the error model - every + // reflection at the floor reads 1/b exactly, however strong it is, and on a merge whose ISa falls + // below the cut NOTHING passes, so every operator is left with no pairs and the point group + // collapses to 1. Measured over the rotation battery: max I/sigma equals 1/b on every merge, and + // the four crystals that lose symmetry are exactly the four whose search merge sits at ISa 3.5-4.0 + // - just above this cut - and drops below 3 when the integration background changes. + // So cap the cut at the merge's own I/sigma quantile: the correlation stage always keeps at least + // its strongest quarter. This is a no-op on any merge where the fixed cut already keeps that many. + constexpr double MIN_PRESENT_FRACTION = 0.25; + double present_cut = opt.present_i_over_sigma; + { + std::vector v; + v.reserve(n); + for (size_t i = 0; i < n; ++i) + if (pass_absence[i]) + v.push_back(IoverSigma[i]); + if (!v.empty()) { + const size_t k = static_cast((1.0 - MIN_PRESENT_FRACTION) * (v.size() - 1)); + std::nth_element(v.begin(), v.begin() + k, v.end()); + present_cut = std::min(present_cut, v[k]); + } + } + + // The correlation stage uses only genuinely-present reflections. Near-zero (systematically + // absent) reflections would otherwise form a second cluster at the origin and fake a high + // correlation for false operators - fatal on centered lattices, where half the reflections + // are extinct. + for (size_t i = 0; i < n; ++i) + pass_cc[i] = pass_absence[i] && IoverSigma[i] >= present_cut && + (opt.min_i_over_sigma <= 0 || IoverSigma[i] >= opt.min_i_over_sigma); + // Resolution-normalised intensity E^2 = I / (shell), from equal-count resolution shells over // the reflections the absence test uses. Lets the absence test judge "present" by intensity // magnitude, not by a possibly under-estimated sigma (see present_e_squared). diff --git a/image_analysis/scale_merge/SearchSpaceGroup.h b/image_analysis/scale_merge/SearchSpaceGroup.h index c0b6fb2e..63ab00c5 100644 --- a/image_analysis/scale_merge/SearchSpaceGroup.h +++ b/image_analysis/scale_merge/SearchSpaceGroup.h @@ -187,6 +187,9 @@ struct SearchSpaceGroupOptions { // Signed I/sigma above which a reflection counts as genuinely present. Signed, so negative // noise is never mistaken for a real reflection. Used both to spot reflections that violate a // wrongly-assumed absence, and to keep the correlation stage from pairing near-zero reflections. + // For the CORRELATION stage this is an upper bound on the cut rather than the cut itself: merged + // sigma is floored at b|I|, so no reflection can read above ISa = 1/b, and a merge whose ISa fell + // below this value would otherwise contribute nothing at all (see SearchSpaceGroup.cpp). double present_i_over_sigma = 3.0; // Extra intensity gate for the systematic-absence test: a reflection also has to reach this