Reindex into the ASU: follow the change of hand, not the hand landed on
ReindexMergedIntoAsu swapped I(+)/I(-) whenever the reindexed index came out on the minus side of its Friedel pair. That is right only for a row that started on the plus side, which is every row when merge_friedel is set - and none of the problem cases. With the mates kept apart, the merge stores the plus mate at +hkl_asu and the minus mate at -hkl_asu (Merge.cpp), while RotationScaleMerge attaches I_plus / I_minus in the plus convention on both rows alike, keying the lookup at plus=true. An alternative-indexing operator is rotation-type, so op(-x) == -op(x) and the two mates always land on opposite sides: testing key.plus alone therefore swaps exactly one row of every pair, whichever way the operator went. The pair ends up self-contradictory, and the mmCIF's pdbx_I_plus / pdbx_I_minus and the F(+)/F(-) columns are swapped on half the rows. Reached by "rugnux -A --model model.pdb" - -A clears MergeFriedel - on a crystal whose indexing ambiguity the model resolves, which is the one path that both keeps the mates apart and reindexes. The swap now fires on was_plus != key.plus. With merge_friedel every row is stored at the ASU representative, so was_plus is always true and this reduces to the previous test; the existing change-of-hand case is unchanged. The test that covered the merge_friedel=false path checked only I, so it passed either way. The new one builds both mates of a pair and asserts they still agree about which intensity is which afterwards - an invariant that holds however the operator moves the hand, so it fails on the old code for either direction. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_016L1qig74oYQzfUJJZbbxFh
This commit is contained in:
@@ -78,7 +78,16 @@ std::vector<MergedReflection> ReindexMergedIntoAsu(const std::vector<MergedRefle
|
||||
for (auto &r : out) {
|
||||
const gemmi::Op::Miller h = op.apply_to_hkl({{r.h, r.k, r.l}});
|
||||
const HKLKey key = key_gen(h[0], h[1], h[2]);
|
||||
if (!key.plus) {
|
||||
// The hands follow a CHANGE of hand, not the hand the new index happens to land on. An
|
||||
// anomalous merge stores the minus mate at -hkl_asu (Merge.cpp), so a row can already be on
|
||||
// the minus side before the operator is applied, while I_plus/I_minus are attached in the
|
||||
// plus convention on both mates alike (RotationScaleMerge.cpp keys the lookup at plus=true).
|
||||
// Since the operator is rotation-type, op(-x) == -op(x), and the minus mate's key.plus is
|
||||
// therefore the negation of the plus mate's - testing key.plus alone swaps exactly the rows
|
||||
// that must not move, and leaves the ones that must. With merge_friedel every row is stored
|
||||
// at the ASU representative, so was_plus is always true and this reduces to the old test.
|
||||
const bool was_plus = key_gen(r.h, r.k, r.l).plus;
|
||||
if (was_plus != key.plus) {
|
||||
std::swap(r.I_plus, r.I_minus);
|
||||
std::swap(r.sigma_plus, r.sigma_minus);
|
||||
std::swap(r.F_plus, r.F_minus);
|
||||
|
||||
@@ -142,3 +142,51 @@ TEST_CASE("Reindex into the ASU: the change of hand swaps the Bijvoet halves", "
|
||||
CHECK(anom[0].l == -asu.l);
|
||||
CHECK(anom[0].I == 100.0f);
|
||||
}
|
||||
|
||||
TEST_CASE("Reindex into the ASU: both mates of an anomalous pair follow the same hand", "[reindex]") {
|
||||
// With the mates kept apart, the merge stores the plus mate at +hkl_asu and the minus mate at
|
||||
// -hkl_asu, while I_plus/I_minus are attached in the plus convention on BOTH rows alike. An
|
||||
// alternative-indexing operator is rotation-type, so op(-x) == -op(x): whatever it does to one
|
||||
// mate's hand it does to the other's, and the two rows must still agree afterwards about which
|
||||
// intensity is which. Deciding the swap from the hand the new index lands on instead of from the
|
||||
// CHANGE of hand moves exactly one of the two, and the pair ends up self-contradictory.
|
||||
const gemmi::SpaceGroup *sg = gemmi::find_spacegroup_by_number(75); // P4, Laue 4/m
|
||||
REQUIRE(sg != nullptr);
|
||||
const HKLKeyGenerator key_gen(/*merge_friedel=*/false, *sg);
|
||||
const HKLKey asu = key_gen(3, 1, 2);
|
||||
REQUIRE(asu.plus);
|
||||
|
||||
// The alternative indexing of a P4 lattice: 4/mmm holds it, 4/m does not.
|
||||
const gemmi::Op op = gemmi::parse_triplet("k,h,-l");
|
||||
|
||||
MergedReflection plus;
|
||||
plus.h = asu.h; plus.k = asu.k; plus.l = asu.l;
|
||||
plus.I = 100.0f; plus.sigma = 1.0f; plus.d = 10.0f;
|
||||
plus.I_plus = 110.0f; plus.sigma_plus = 2.0f; plus.I_minus = 90.0f; plus.sigma_minus = 3.0f;
|
||||
plus.F_plus = 10.5f; plus.F_minus = 9.5f;
|
||||
|
||||
MergedReflection minus = plus; // same anomalous split, stored at -hkl_asu
|
||||
minus.h = -asu.h; minus.k = -asu.k; minus.l = -asu.l;
|
||||
|
||||
const auto out = ReindexMergedIntoAsu({plus, minus}, op, 75, /*merge_friedel=*/false);
|
||||
REQUIRE(out.size() == 2);
|
||||
|
||||
// The two rows are still one pair: same Laue label up to the Friedel sign, opposite hands.
|
||||
CHECK(out[0].h == -out[1].h);
|
||||
CHECK(out[0].k == -out[1].k);
|
||||
CHECK(out[0].l == -out[1].l);
|
||||
|
||||
// ... and they still tell the same story about which hand carries which intensity.
|
||||
CHECK(out[0].I_plus == out[1].I_plus);
|
||||
CHECK(out[0].I_minus == out[1].I_minus);
|
||||
CHECK(out[0].sigma_plus == out[1].sigma_plus);
|
||||
CHECK(out[0].sigma_minus == out[1].sigma_minus);
|
||||
CHECK(out[0].F_plus == out[1].F_plus);
|
||||
CHECK(out[0].F_minus == out[1].F_minus);
|
||||
|
||||
// And the swap did have to happen: "k,h,-l" takes (3,1,2) to (1,3,-2), whose Laue class has its
|
||||
// ASU representative at l > 0, so the row that was the plus mate of its pair is the minus mate of
|
||||
// the new one. The intensity now standing at the plus index is the one that was at -hkl before.
|
||||
CHECK(out[0].I_plus == 90.0f);
|
||||
CHECK(out[0].I_minus == 110.0f);
|
||||
}
|
||||
|
||||
Reference in New Issue
Block a user