Merge branch 'd-5ebi' into rc172-2

image_analysis/scale_merge/SearchSpaceGroup.cpp

@@ -1377,6 +1377,8 @@ SearchSpaceGroupResult SearchSpaceGroup(
                     // reason as h_ratio: to report the adopted candidate's whether or not it decided.
                     double r_added = std::numeric_limits<double>::quiet_NaN();
                     double r_over_floor = std::numeric_limits<double>::quiet_NaN();
+                    // The candidate's OWN operators, worst R over best R (see max_operator_r_spread).
+                    double r_spread = std::numeric_limits<double>::quiet_NaN();
                     // Report-only, for the finalist ledger: the subgroup b the tests below compared
                     // this candidate against, whether it survived every consistency test, and - when
                     // it was the highest refusal - why it did not. The selection loop computes all
@@ -1629,6 +1631,23 @@ SearchSpaceGroupResult SearchSpaceGroup(
             if (n_new > 0 && pairs_new >= opt.min_pairs_for_h)
                 r_added = r_new / n_new;
         }
+        // Homogeneity of the candidate's OWN operators: the worst-agreeing over the best-agreeing,
+        // over all of them, parent and added alike (see max_operator_r_spread).
+        {
+            double r_lo = std::numeric_limits<double>::infinity(), r_hi = 0.0;
+            int n_op = 0;
+            for (const auto &rot : c.pg->rotations) {
+                if (rot.rot == gemmi::Op::identity().rot)
+                    continue;
+                const auto &os = operator_score(rot);
+                if (os.n_pairs < opt.min_pairs_per_operator || !std::isfinite(os.r_stat)
+                    || os.r_stat <= 0.0)
+                    continue;
+                r_lo = std::min(r_lo, os.r_stat); r_hi = std::max(r_hi, os.r_stat); ++n_op;
+            }
+            if (n_op >= 2)
+                c.r_spread = r_hi / r_lo;
+        }
         c.parent_b_used = parent_b;
         c.r_added = r_added;
         c.r_over_floor = (std::isfinite(r_added) && std::isfinite(result.merge_r_floor)
@@ -1686,16 +1705,28 @@ SearchSpaceGroupResult SearchSpaceGroup(
         if (r_refused)
             consistent = false;
 
+        // Operator-homogeneity test. A necessary condition like H and the contrast, and the one that
+        // catches what a MEAN over the added operators cannot: a candidate whose added set mixes real
+        // operators with false ones. The gates above read summaries - mean H, mean R - so a promotion
+        // that adds two operators at the noise floor together with four several times worse reads a
+        // mean between them and passes, while the group it names is a mixture that cannot be a
+        // symmetry of one crystal. Asking instead that the candidate's operators agree with EACH OTHER
+        // sees the mixture directly (see max_operator_r_spread).
+        const bool spread_refused = std::isfinite(c.r_spread) && c.r_spread > opt.max_operator_r_spread;
+        if (spread_refused)
+            consistent = false;
+
         if (!consistent) {
             // Same three sentences the highest refusal gets below, but kept per candidate for the
             // ledger: short, because the long form is written once for the group the user is told about.
-            c.why = h_refused ? "H ratio" : r_refused ? "added-operator R" : "merge chi^2";
+            c.why = h_refused ? "H ratio" : r_refused ? "added-operator R"
+                  : spread_refused ? "operator spread" : "merge chi^2";
             // Record the highest-order refusal so the caller can say WHY it is processing lower.
             if (c.order > refused_order && c.pg->representative) {
                 refused_order = c.order;
                 refused_pg_hm = c.pg->representative->point_group_hm();
                 refused_pg_rep = c.pg->representative;
-                refused_twin_gate = h_refused || r_refused;
+                refused_twin_gate = h_refused || r_refused || spread_refused;
                 refused_fraction = h_refused ? std::max(0.0, 0.5 - h_added)
                                              : std::numeric_limits<double>::quiet_NaN();
                 if (h_refused)
@@ -1727,6 +1758,13 @@ SearchSpaceGroupResult SearchSpaceGroup(
                               + FormatDouble(opt.max_operator_r_over_floor, 1)
                               + "x) - at full resolution the added operator relates reflections that "
                                 "disagree far beyond measurement error, i.e. it is not a real symmetry";
+                else if (spread_refused)
+                    refused_why = "its own operators disagree among themselves - the worst-agreeing "
+                                  "reads " + FormatDouble(c.r_spread, 1)
+                                + "x the intensity-weighted R of the best-agreeing one in the SAME "
+                                  "group (bound " + FormatDouble(opt.max_operator_r_spread, 1)
+                                + "x) - so the group pools operators the data support with operators "
+                                  "they do not, and cannot be a symmetry of one crystal";
                 else if (std::isfinite(c.chi2) && std::isfinite(chi2_ref))
                     refused_why = "merge chi^2 is " + FormatDouble(c.chi2 / chi2_ref, 2)
                                 + "x the subgroup's (bound " + FormatDouble(opt.max_merge_chi2_ratio, 2) + ")";
@@ -1822,6 +1860,7 @@ SearchSpaceGroupResult SearchSpaceGroup(
                                / (result.random_pairing_r - global_best_r);
         }
         e.r_over_floor = c.r_over_floor;
+        e.r_spread = c.r_spread;
         e.adopted = (c.pg == best_pg);
         e.eligible = c.eligible;
         e.refused_reason = c.why;
@@ -2589,14 +2628,15 @@ std::string FinalistLedgerToText(const SearchSpaceGroupResult& result) {
         os << "  Unrelated reflections (shell-matched random pairing) R "
            << FormatDouble(result.random_pairing_r, 4)
            << " - the far end of the contrast scale.\n";
-    os << "\n  PG      ord  minCC   R_add contrast R/best R/floor  |  H     chi2/best   b     b/par  verdict\n"
-          "                 <---- calibrated evidence ---->        |  <--- diagnostic only, NOT evidence --->\n";
+    os << "\n  PG      ord  minCC   R_add contrast spread R/best R/floor  |  H     chi2/best   b     b/par  verdict\n"
+          "                 <------- calibrated evidence ------->        |  <--- diagnostic only, NOT evidence --->\n";
     for (const auto& e : result.point_group_ledger) {
         os << "  " << std::left << std::setw(8) << e.point_group_hm << std::right << std::setw(3)
            << e.order
            << std::setw(8) << LedgerCell(e.min_class_cc, 3)
            << std::setw(8) << LedgerCell(e.r_added, 4)
            << std::setw(9) << LedgerCell(e.r_contrast, 2)
+           << std::setw(7) << LedgerCell(e.r_spread, 2)
            << std::setw(7) << LedgerCell(e.r_over_best, 2)
            << std::setw(8) << LedgerCell(e.r_over_floor, 2)
            << "  |" << std::setw(6) << LedgerCell(e.h_ratio, 2)

image_analysis/scale_merge/SearchSpaceGroup.h

@@ -471,6 +471,26 @@ struct SearchSpaceGroupOptions {
     //    the floor is a noisier reference than a confirmed operator (it can also be inflated by an
     //    indexing ambiguity that mixes hands), and the safe direction here is to stay in P1.
     double max_operator_r_over_floor = 3.5;
+    //
+    //  - max_operator_r_spread: the two tests above read a MEAN over the operators a promotion adds,
+    //    so they are blind to an added set that MIXES real operators with false ones - the mean sits
+    //    between the two populations and passes. That is what a pseudo-symmetric metric produces: a
+    //    monoclinic crystal whose lattice is exactly hexagonal (a = c, beta = 120) carries two genuine
+    //    horizontal 2-folds among the six a 622 promotion adds, and their R at the noise floor buys
+    //    the other four a mean well inside the contrast bound. The mixture is invisible to any summary
+    //    of the added set and obvious the moment the candidate's operators are compared with EACH
+    //    OTHER: a point group is the claim that all of its operators are symmetries of one crystal, so
+    //    they must agree with the merge about equally well.
+    //
+    //    Measured over the 129 rotation datasets of one battery run that adopted a group of order 3 or
+    //    more, as worst-agreeing over best-agreeing operator WITHIN the adopted group: every correct
+    //    adoption sits at 1.0-2.64, the highest being the weak cubic crystal whose near-spindle 3-fold
+    //    is unusually clean (the same crystal that motivated the contrast above), then 2.30 for a
+    //    fine-sliced tetragonal one. The one over-called point group on that run reads 4.87. 3.2 sits
+    //    in that gap with 0.56 of headroom below and 1.67 above - deliberately well clear of the
+    //    2.0-2.6 spread a genuine group reaches, which is exactly where the old ratio-to-best gate's
+    //    bound of 2.0 sat and why that gate had to go.
+    double max_operator_r_spread = 3.2;
 
     // 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
@@ -680,6 +700,7 @@ struct PointGroupLedgerEntry {
     double r_over_best = std::numeric_limits<double>::quiet_NaN();     // ... over the best-agreeing operator anywhere
     double r_contrast = std::numeric_limits<double>::quiet_NaN();      // ... 0 = unrelated reflections, 1 = that best operator
     double r_over_floor = std::numeric_limits<double>::quiet_NaN();    // ... over the merge's random-noise R floor
+    double r_spread = std::numeric_limits<double>::quiet_NaN();        // this group's own worst operator R over its best
     bool adopted = false;         // this is the group the search took
     bool eligible = false;        // it passed every consistency test (so it COULD have been taken)
     std::string refused_reason;   // why it did not, when it was the highest refusal