rugnux: fix de-novo space-group centering and weak-data high symmetry
De-novo space-group search: - Determine lattice centering from the strength of the systematically-absent class (its mean I/sigma vs the present class) instead of a per-reflection violation count. The count gate was brittle on noisy / obverse-reverse- twinned data, where enough genuinely-absent reflections cross I/sigma>3 to trip the 10% bound although the class is several-fold weaker, and wrongly kept a genuine R-centered lattice primitive. False centerings, whose absent class is as strong as the present one, stay rejected. Screws keep the count test (too few axial reflections to average). - Recover a genuine high-symmetry point group when the merged error model is badly miscalibrated (weak, low-resolution data whose sigmas are far too small). The fixed-sigma chi^2 ratio then grows with point-group order for real symmetry too and wrongly rejects it (a true cubic group reached a chi^2 ratio ~14); the self-normalising systematic-b test stays valid, so once chi2_ref shows the ratio is untrustworthy the promotion is confirmed on the b-test alone. The balloon veto still rejects a twin. This is the gap left when the earlier log10(chi2_ref) widening was removed for the b-veto. Also print the adopted space group and unit cell as a clean one-line summary at the end of the run (de-novo or user-fixed -S), not only in the mmCIF. Full rotation-test battery: exactly the two intended crystals change (both to the correct answer, with lower R-meas and higher ISa), every other crystal unchanged; space-group match 21->22/25. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
This commit is contained in:
@@ -5,6 +5,9 @@ This is an UNSTABLE release. It includes many experimental features, as well as
|
||||
|
||||
* rugnux: Add a dataset-wide **Wilson B-factor** estimate to the merged output (mmCIF `_reflns.B_iso_Wilson_estimate`, the printed merge statistics and the log) - the analogue of XDS's Wilson-line B, which was not exported. Fitted from the shell-averaged intensity vs 1/d^2 over the meaningful resolution range; diagnostic only, not fed back into scaling. The per-image Wilson B (viewer plot) now emits NaN for an implausible fit instead of the spurious hundreds-of-A^2 value a bad frame used to produce.
|
||||
* rugnux: Veto a merohedral-twin over-promotion in the de-novo space-group search. A partial twin whose merged chi^2 looks self-consistent could be promoted to the holohedral group even though merging its non-equivalent reflections balloons the error-model b; a chi^2-passing promotion whose systematic b balloons past a calibrated bound is now kept in its true lower symmetry (an R3 case previously merged as R32; battery space-group match 20->21/25, no change to any genuine high-symmetry crystal).
|
||||
* rugnux: De-novo space-group search - decide lattice **centering** from the strength of the systematically-absent class instead of a per-reflection violation count. A real centering cancels structure factors, so its absent class is systematically weak (its mean I/sigma sits well below the present class); the previous test - how many absent reflections individually cleared I/sigma>3 - was brittle on noisy or obverse/reverse-twinned data, where enough genuinely-absent reflections cross the threshold to trip the bound although the class is several-fold weaker. A rhombohedral (R) lattice whose absent class ran to ~13% individual violations was wrongly kept primitive (P) and is now correctly centered (with lower R-meas and higher ISa, matching the reference); the false-centering cases stay rejected because their absent class is as strong as the present one. Screw/glide axes keep the count test - their predicted-absent class is a handful of axial reflections, too few to average. Battery space-group match 21->22/25.
|
||||
* rugnux: De-novo space-group search - recover a genuine high-symmetry point group on weak, low-resolution data whose merged error model is badly miscalibrated. When even the best subgroup's reduced chi^2 is far above 1 (merged sigmas badly under-estimated), the fixed-sigma chi^2 ratio grows with point-group order for real symmetry too and wrongly rejects it (a true cubic group reached a chi^2 ratio of ~14). The systematic-b test re-fits its own error and stays valid under a broken sigma model, so once the chi^2 ratio is shown untrustworthy the promotion is confirmed on the b-test alone - a genuine step's b barely moves, whereas a twin's balloons and is still vetoed. Restores a weak low-resolution cubic (F432) crystal that had fallen to a tetragonal subgroup, with no change to any other crystal in the battery. (This is the gap left when the earlier log10(chi2_ref) chi^2 widening was removed in favour of the b-veto; the b-test now covers it without re-admitting the twin.)
|
||||
* rugnux: Print the adopted **space group and unit cell** as a one-line summary at the end of the run (`Space group: <name> (No. N)` alongside `Unit cell:`), for both a de-novo determination and a user-fixed `-S` space group - previously the determined space group appeared only inside the de-novo search diagnostics or in the mmCIF.
|
||||
* rugnux: Score the radiation-damage **decay** cross-validation on the same sigma-independent (R-meas-like) metric as the absorption/modulation surfaces, so a spurious decay slope can no longer pass the gate by reshaping the sigmas rather than tightening the intensities.
|
||||
* rugnux: Fix de-novo rotation indexing occasionally committing a spurious axis-multiple supercell (a true axis reported at 2x/3x its length), which collapsed the whole dataset to P1. Two-pass indexing scored its first-pass schemes by how many validation frames each indexes, but an axis multiple indexes every frame its true cell does, so both schemes saturate and the tie fell to whichever ran first; when that scheme's full-rotation FFT resolved a true axis only as its harmonic, the multiplied cell was committed. A cross-scheme tie-break now prefers the smaller cell whenever the two schemes' cell volumes are related by a near-integer factor (the signature of a spurious axis multiple, independent of scheme order), while a genuine centering relation (a non-integer volume ratio) is left untouched; the same axis-doubling is additionally unmasked when selecting among a single scheme's own candidate cells, where a sub-cell that still indexes essentially as many spots overrides a cell more than 1.5x larger (its own primitive).
|
||||
* rugnux: Widen the refined-cell angle bounds to [30, 150] degrees for both rotation candidate refinement and per-frame stills cell refinement, so a genuinely oblique lattice is no longer clamped towards 90 degrees, and check a refined candidate's angles against the reference cell in `Refine()` when a reference is supplied.
|
||||
|
||||
@@ -482,9 +482,18 @@ SearchSpaceGroupResult SearchSpaceGroup(
|
||||
parent_b = s.b_extra;
|
||||
}
|
||||
}
|
||||
if (!consistent && std::isfinite(c.chi2) && std::isfinite(chi2_ref)
|
||||
&& c.chi2 <= chi2_ref * opt.max_merge_chi2_rescue
|
||||
&& parent_b > 1e-4 && c.b_extra <= parent_b * opt.max_systematic_b_ratio)
|
||||
// The chi^2 ratio is only trustworthy when the error model is calibrated. When even the best
|
||||
// subgroup's reduced chi^2 (chi2_ref) is far above 1 - weak, low-resolution data whose merged
|
||||
// sigmas are badly under-estimated - the ratio grows with point-group order for genuine high
|
||||
// symmetry too and wrongly rejects it (a true weak F432 reaches ratio ~14). The systematic-b test
|
||||
// re-fits its own error, so it stays valid under a broken sigma model: a genuine step's b barely
|
||||
// moves (b-ratio ~1) while a twin's balloons. So once chi2_ref shows the error model is unreliable,
|
||||
// a promotion is rescued on the b-test alone (subject to the balloon veto below); otherwise the
|
||||
// rescue is confined to the narrow chi^2 band just past the ratio bound.
|
||||
const bool miscalibrated = std::isfinite(chi2_ref) && chi2_ref > opt.chi2_ref_reliable;
|
||||
if (!consistent && parent_b > 1e-4 && c.b_extra <= parent_b * opt.max_systematic_b_ratio
|
||||
&& (miscalibrated || (std::isfinite(c.chi2) && std::isfinite(chi2_ref)
|
||||
&& c.chi2 <= chi2_ref * opt.max_merge_chi2_rescue)))
|
||||
consistent = true;
|
||||
// Veto a chi^2-passing promotion whose b clearly ballooned (above the largest genuine step, below a
|
||||
// twin); a genuine but imperfectly-scaled high-symmetry merge stays under the bound and is untouched.
|
||||
@@ -531,7 +540,8 @@ SearchSpaceGroupResult SearchSpaceGroup(
|
||||
// Judge centering and screw/glide absences on separate reflection sets. Lumping them lets
|
||||
// a large, correct centering-absent set hide a few strong screw violations and over-claim
|
||||
// screw axes (e.g. I4_132 on I432 data).
|
||||
int centering_absent = 0, centering_violations = 0;
|
||||
int centering_absent = 0;
|
||||
double centering_absent_sum = 0;
|
||||
int screw_absent = 0, screw_violations = 0;
|
||||
|
||||
for (size_t i = 0; i < n; ++i) {
|
||||
@@ -549,7 +559,8 @@ SearchSpaceGroupResult SearchSpaceGroup(
|
||||
s.absent_observed += 1;
|
||||
absent_sum += IoverSigma[i];
|
||||
centering_absent += 1;
|
||||
if (present) { s.absent_violations += 1; centering_violations += 1; }
|
||||
centering_absent_sum += IoverSigma[i];
|
||||
if (present) s.absent_violations += 1;
|
||||
} else if (gops.is_systematically_absent(hkl)) {
|
||||
s.absent_observed += 1;
|
||||
absent_sum += IoverSigma[i];
|
||||
@@ -566,8 +577,20 @@ SearchSpaceGroupResult SearchSpaceGroup(
|
||||
if (present_n > 0)
|
||||
s.present_mean_i_over_sigma = present_sum / present_n;
|
||||
|
||||
// Centering is judged by class STRENGTH, not a per-reflection violation count. A real centering
|
||||
// cancels structure factors, so its absent class is systematically weak - its mean signed
|
||||
// I/sigma sits well below the present class - regardless of noise or obverse/reverse twinning;
|
||||
// a false centering leaves the "absent" class as strong as the present one (mean ratio ~1). The
|
||||
// count-of-strong-violations gate is brittle on noisy/twinned data, where enough genuinely-absent
|
||||
// reflections randomly clear I/sigma>3 to trip the 10% bound though the class is 3-4x weaker (a
|
||||
// true R3 at 13.5% violations, absent 1.7 vs present 6.0). The mean is well-determined here
|
||||
// because a centering-absent class holds a third-to-half of all reflections. Screws keep the
|
||||
// count gate: their predicted-absent class is a handful of axial reflections, too few to average.
|
||||
const double present_mean = present_n > 0 ? present_sum / present_n : 0.0;
|
||||
const double centering_absent_mean =
|
||||
centering_absent > 0 ? centering_absent_sum / centering_absent : 0.0;
|
||||
const bool centering_ok = centering_absent == 0 ||
|
||||
centering_violations <= opt.max_absent_violation_fraction * centering_absent;
|
||||
(present_n > 0 && centering_absent_mean <= opt.max_absent_present_ratio * present_mean);
|
||||
const bool screw_ok = screw_absent == 0 ||
|
||||
screw_violations <= opt.max_absent_violation_fraction * screw_absent;
|
||||
s.consistent = centering_ok && screw_ok;
|
||||
|
||||
@@ -117,6 +117,14 @@ struct SearchSpaceGroupOptions {
|
||||
// keeps a clearly-ballooned promotion down, never promotes.
|
||||
double max_systematic_b_veto = 2.0;
|
||||
|
||||
// Above this reduced chi^2 for the best subgroup (chi2_ref), the merged error model is treated as
|
||||
// badly miscalibrated (weak, low-resolution data whose sigmas are far too small): the fixed-sigma
|
||||
// chi^2 ratio then grows with point-group order for genuine high symmetry too and can no longer
|
||||
// arbitrate, so a promotion is confirmed on the systematic-b test alone (which re-fits its own error
|
||||
// and stays valid). A well-calibrated merge sits near 1; 3.0 (sigmas ~1.7x too small) marks the point
|
||||
// where the ratio stops being trustworthy. The balloon veto still guards against a twin.
|
||||
double chi2_ref_reliable = 3.0;
|
||||
|
||||
// --- Stage B: space group (screw axes / centering) ---
|
||||
bool determine_space_group = true; // false: stop at the symmorphic representative
|
||||
|
||||
@@ -132,10 +140,19 @@ struct SearchSpaceGroupOptions {
|
||||
// (e.g. the monoclinic 2_1: 0k0-odd at ~1% of 0k0-even). 0 disables the gate (I/sigma only).
|
||||
double present_e_squared = 0.3;
|
||||
|
||||
// A candidate's absence conditions are accepted when at most this fraction of the reflections
|
||||
// it predicts absent are in fact strongly present.
|
||||
// A candidate's SCREW/glide absence conditions are accepted when at most this fraction of the
|
||||
// reflections it predicts absent are in fact strongly present.
|
||||
double max_absent_violation_fraction = 0.10;
|
||||
|
||||
// A CENTERING is accepted when its systematically-absent class is this much weaker than the present
|
||||
// class - mean signed I/sigma of the centering-absent reflections <= this fraction of the present
|
||||
// mean. A real centering cancels structure factors so its absent class sits near zero (ratio ~0-0.3
|
||||
// across the test battery) even on noisy or obverse/reverse-twinned data; a false centering leaves it
|
||||
// as strong as the present class (ratio ~1.0). 0.5 separates the two with wide margin. This strength
|
||||
// test replaces a per-reflection violation-count gate for centering, which was brittle when noise
|
||||
// pushed genuinely-absent reflections over I/sigma>3 (a true R3 was lost at 13.5% violations).
|
||||
double max_absent_present_ratio = 0.5;
|
||||
|
||||
// Need at least this many observed reflections in the predicted-absent class before a
|
||||
// screw/centering is claimed (guards against deciding from a handful of reflections).
|
||||
int min_absent_observed = 8;
|
||||
|
||||
@@ -1517,6 +1517,14 @@ int main(int argc, char **argv) {
|
||||
std::cout << fmt::format("Total throughput:{:.2f} MB/s", result.throughput_MBs) << std::endl;
|
||||
if (result.indexing_rate.has_value())
|
||||
std::cout << fmt::format("Indexing rate: {:.2f}%", result.indexing_rate.value() * 100.0) << std::endl;
|
||||
// Final one-line summary of the adopted crystal (whether de-novo determined or fixed with -S),
|
||||
// so it is not buried in the space-group-search block (which is de-novo only) or only in the mmCIF.
|
||||
if (result.space_group_number.has_value()) {
|
||||
const gemmi::SpaceGroup *sg = gemmi::find_spacegroup_by_number(
|
||||
static_cast<int>(result.space_group_number.value()));
|
||||
std::cout << fmt::format("Space group: {} (No. {})", sg ? sg->short_name() : "?",
|
||||
result.space_group_number.value()) << std::endl;
|
||||
}
|
||||
if (result.consensus_cell.has_value()) {
|
||||
const auto &c = result.consensus_cell.value();
|
||||
std::cout << fmt::format("Unit cell: a={:.2f} b={:.2f} c={:.2f} alpha={:.2f} beta={:.2f} gamma={:.2f}",
|
||||
|
||||
Reference in New Issue
Block a user