diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index dfcc7b37..b3bb229d 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -3,6 +3,7 @@ ### 1.0.0-rc.165 This is an UNSTABLE release. It includes many experimental features, as well as many AI generated fixes. We recommend using rc.152 for production use. +* `MEAN_ATOM_DENSITY_SIGMA`, the mean 2Fo-Fc density at the model's atom centres that `rugnux --model` reports, is read from the map by cubic rather than linear interpolation and comes out around a tenth higher as a result. The map, the R-factors and every other number are unchanged; only this reading of the map moves, so it is no longer comparable with the figure earlier versions printed. * `rugnux --model` writes an anomalous difference map (`_anom.ccp4`) when the merge kept the Bijvoet split, and names the ten model atoms it peaks highest on - `ANOMALOUS_SITE_01`..`_10` in the results report, and one line in the log - so a run says which atoms carry the anomalous signal and how many sigma each one is. * The rugnux results report records how the run was invoked, what it cost and what it ran on: `COMMAND_LINE=` is the command line as one shell-ready line, `WALL_TIME=` the whole invocation in seconds, and `GPU_COUNT=` / `GPU=` how many GPUs were visible and which ones. The total wall time is also printed on stdout, next to the processing time it is slightly larger than. * The rugnux results report has a `10. MODEL VALIDATION` section when `--model` was given: R-work / R-free with their reflection counts, the bulk-solvent and scale parameters, the mean 2Fo-Fc density at the atom centres, the reindexing operators and where the maps were written - or `MODEL_VALIDATION= NOT_PERFORMED` and the reason, when the model could not be used. `REPORT_VERSION` is 4; `WARNINGS` moves from section 10 to 11 and no existing key changed. diff --git a/rugnux/ModelValidation.cpp b/rugnux/ModelValidation.cpp index 1072650d..3387c94d 100644 --- a/rugnux/ModelValidation.cpp +++ b/rugnux/ModelValidation.cpp @@ -49,6 +49,11 @@ double write_ccp4(const gemmi::Grid &map, const std::string &path) { return ccp4.hstats.rms; } +// Cubic, not the default linear, for reading a map at a point. The maps are sampled every d_min/3, +// and a peak that sharp read by trilinear interpolation comes out up to a quarter low - unevenly +// enough to reorder the anomalous sites. +constexpr int MAP_INTERPOLATION_ORDER = 3; + // How many anomalous sites the report names. The strongest few are what says whether the anomalous // signal is there and what carries it; a full site list is what the map file is for. constexpr size_t MAX_ANOMALOUS_SITES = 10; @@ -285,7 +290,7 @@ ModelValidationResult ValidateAgainstModel(const std::vector & for (gemmi::Model &m : st.models) for (gemmi::Chain &ch : m.chains) for (gemmi::Residue &r : ch.residues) - for (gemmi::Atom &a : r.atoms) { s += grid2fofc.interpolate_value(a.pos); ++n; } + for (gemmi::Atom &a : r.atoms) { s += grid2fofc.interpolate_value(a.pos, MAP_INTERPOLATION_ORDER); ++n; } result.mean_atom_density_sigma = (n > 0 && rms2 > 0) ? (s / n) / rms2 : 0; } @@ -347,12 +352,10 @@ ModelValidationResult ValidateAgainstModel(const std::vector & for (gemmi::Atom &a : r.atoms) { if (a.is_hydrogen()) // hydrogen scatters no anomalous signal continue; - // Cubic, not the default linear: the map is sampled every d_min/3, and a - // peak that sharp read by trilinear interpolation comes out up to a quarter - // low - unevenly enough to reorder the sites. sites.push_back({fmt::format("{} {} {}{}", a.name, r.name, ch.name, r.seqid.str()), - rms > 0 ? grid.interpolate_value(a.pos, 3) / rms : 0.0}); + rms > 0 ? grid.interpolate_value(a.pos, MAP_INTERPOLATION_ORDER) / rms + : 0.0}); } std::sort(sites.begin(), sites.end(), [](const auto &x, const auto &y) { return x.sigma > y.sigma; });