Solvent build: write each new water beside the model, not a symmetry mate

A difference peak is taken once per symmetry-equivalent set (the map's
asymmetric unit) and was to be moved to its copy nearest the model's own
atoms (NearestCopy) before a water is built there. That step took the
nearest atom's NeighborSearch mark, whose position is already the
symmetry image, and applied the image's operator to it a second time:
outside the identity image the water was written beside a symmetry mate
or a unit cell away (by a lattice translation for a two-fold, beside
another mate for a three-, four- or six-fold).

The copy is now taken against the atom as the model has it: the
symmetry operator and lattice translation that put the peak nearest to
that atom (UnitCell::find_nearest_image / find_nearest_pbc_position). The
water's environment, and so the refinement, is the same; only the
coordinates written move. New test: a peak beside a four-fold mate in
P4 is written 2.8 A from the model's water (the old step wrote it a
lattice vector away).

5reo (open arm, C2), --refine-solvent, share of the added waters whose
nearest non-water, non-hydrogen model atom over all copies is in the
written copy: from the model without waters 67/155 (43%) before,
154/155 after; from the deposited model 27/43 before, 42/43 after
(the one left a near tie moved by refinement). R-work/R-free
0.2046/0.2335 -> 0.2046/0.2334 and 0.1979/0.2172 unchanged; 156 and 43
waters added either way.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01UqXYdL4qBGYRiSXSGXfCLu

docs/CPU_DATA_ANALYSIS_REFINEMENT.md

@@ -476,7 +476,9 @@ higher, one 0.0016 lower). `--refine-score-only` builds nothing, and the build a
   restrains the model's own (§15.5).
 
 - *Candidates* are the local maxima of mFo-DFc from 3 σ (phenix.refine's peak cutoff), each symmetry-equivalent
-  set once, moved to the copy nearest the model. A peak closer than 1.2 Å to a model atom is that atom's own
+  set once, moved to the copy nearest the model - the symmetry operator and lattice translation that put it
+  closest to an atom of the model as written, so a water built there is written beside the model, not beside a
+  symmetry mate or in another unit cell. A peak closer than 1.2 Å to a model atom is that atom's own
   residual and is passed over.
 - *Not solvent* - reported in `<prefix>_solvent.json` with the reason, never filled: a peak 1.2-1.9 Å from a model
   atom (a missing or displaced atom); one beside a model atom (within 2.6 Å) that sits in mFo-DFc of -3 σ or

image_analysis/structure_refinement/SolventBuild.cpp

@@ -132,18 +132,15 @@ std::vector<Peak> DifferencePeaks(const gemmi::Grid<float> &g, double rms, doubl
     return out;
 }
 
-// The copy of `p` nearest to the model's own atoms (not their symmetry mates); false where no atom is within
-// `reach`.
+// The copy of `p` - the symmetry operator and lattice translation - nearest to the model's own atoms (not their
+// symmetry mates); false where no atom is within `reach`. The mark's position is already the atom's image, so the
+// copy is taken by the atom as the model has it: p's copy nearest to that atom.
 bool NearestCopy(gemmi::NeighborSearch &ns, const gemmi::UnitCell &cell, double reach, gemmi::Position &p) {
     gemmi::NeighborSearch::Mark *m = ns.find_nearest_atom(p, reach);
     if (m == nullptr)
         return false;
-    const gemmi::FTransform t = ns.get_image_transformation(m->image_idx);
-    const gemmi::Fractional fc = cell.fractionalize(p);
-    const gemmi::Fractional fi = t.apply(cell.fractionalize(m->pos));
-    const gemmi::Fractional shifted(fc.x - std::round(fc.x - fi.x), fc.y - std::round(fc.y - fi.y),
-                                    fc.z - std::round(fc.z - fi.z));
-    p = cell.orthogonalize(gemmi::FTransform(t.inverse()).apply(shifted));
+    const gemmi::Position atom = m->to_cra(*ns.model).atom->pos;
+    p = cell.find_nearest_pbc_position(atom, p, cell.find_nearest_image(atom, p, gemmi::Asu::Any).sym_idx);
     return true;
 }
 

tests/SolventBuildTest.cpp

@@ -101,6 +101,25 @@ TEST_CASE("SolventBuild_peak_with_and_without_partner", "[Solvent]") {
     CHECK(r.edits.at(0).site == "HOH W1");
 }
 
+// A peak beside a symmetry mate of the model's water, in another unit cell, is written beside the water itself: the
+// copy of the peak - symmetry operator and lattice translation - nearest to the model's own atoms.
+TEST_CASE("SolventBuild_peak_written_beside_the_model", "[Solvent]") {
+    gemmi::Structure st = box();
+    st.spacegroup_hm = "P 4";
+    st.setup_cell_images();
+    add_atom(st, "HOH", "O", gemmi::El::O, {6, 4, 6});
+    // the water's mate by -y,x,z is at 20 6 6; the peak 2.8 A from it, only there (the map drawn without symmetry)
+    const gemmi::Position peak(20, 6, 8.8);
+    const gemmi::Grid<float> diff = map_with(gemmi::UnitCell(EDGE, EDGE, EDGE, 90, 90, 90), {peak}, 6.0f);
+    const gemmi::Grid<float> map = map_with(st.cell, {}, 0.0f, 1.5f);
+    const SolventMaps maps{&map, 1.0, &diff, 1.0, nullptr, 0};
+
+    const SolventRound r = PickSolvent(st, maps, {}, {}, {}, 300, 1);
+    REQUIRE(r.added == 1);
+    const SolventEdit &e = r.edits.at(0);
+    CHECK(e.pos.dist(gemmi::Position(6, 4, 6)) == Catch::Approx(2.8).margin(0.25));
+}
+
 // A sulfate-like cluster - a peak with four more at 1.5 A on a tetrahedron - beside a water is no water: the centre
 // is reported as a tetrahedral group, its corners as atoms of one molecule; nothing is added.
 TEST_CASE("SolventBuild_tetrahedral_cluster_refused", "[Solvent]") {