Merge twin-gate into rc181: twin decision against a per-pair null (rule A) and taken at 2-3 SE where the L-test reads twinned (rule B)
Build Packages / Create release (push) Successful in 7s
Build Packages / Build documentation (push) Successful in 47s
Build Packages / build:portable:macos-arm64:nocuda (push) Successful in 5m42s
Build Packages / build:rugnux:linux-aarch64:cuda (push) Successful in 11m32s
Build Packages / build:portable:linux-x86_64:nocuda (push) Successful in 23m19s
Build Packages / build:portable:windows-x86_64:nocuda (push) Failing after 23m47s
Build Packages / build:portable:linux-x86_64:cuda (push) Successful in 26m16s
Build Packages / build:jfjoch:ubuntu2204:nocuda (push) Successful in 26m4s
Build Packages / build:portable:windows-x86_64:cuda (push) Failing after 28m5s
Build Packages / build:jfjoch:rocky8:nocuda (push) Successful in 30m53s
Build Packages / build:jfjoch:rocky9:nocuda (push) Successful in 31m0s
Build Packages / HDF5 consumer tests (DIALS, XDS) (push) Successful in 33m41s
Build Packages / build:jfjoch:ubuntu2404:nocuda (push) Successful in 24m30s
Build Packages / Generate python client (push) Successful in 30s
Build Packages / build:jfjoch:rocky8:cuda-sls9 (push) Successful in 27m51s
Build Packages / build:jfjoch:ubuntu2204:cuda (push) Successful in 26m48s
Build Packages / build:jfjoch:rocky8:cuda (push) Successful in 31m13s
Build Packages / build:jfjoch:rocky9:cuda-sls9 (push) Successful in 33m44s
Build Packages / build:jfjoch:ubuntu2404:cuda (push) Successful in 29m20s
Build Packages / build:jfjoch:rocky9:cuda (push) Successful in 33m9s
Build Packages / Unit tests (push) Successful in 1h36m58s
Build Packages / Validation Battery (push) Successful in 51m10s

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

# Conflicts:
#	rugnux/RugnuxScaleMerge.cpp
This commit is contained in:
leonarski_f committed 2026-10-28 17:53:32 +01:00
commit b13f1bfd2d
8 files changed
+145 -59

No files matched your search

+21 -9
View File
@@ -539,16 +539,28 @@ is linear in the two intensities, so an unrelated partner takes a coefficient of
least-squares fit of the twinned amplitudes instead took a fraction of 0.16 for unrelated partners of one model at
R-free 0.31, the average of two predictions being a better prediction than either.
**The decision.** The twinned refinement is kept where its fraction is above 0, its R-free is below the
untwinned refinement's by more than 0.001, and the law beats its own **null** on the free reflections by at least
three standard errors - the null being the same twinned model at the same fraction with every reflection's mate
replaced by an unrelated reflection of the same resolution: what the averaging of two intensities buys whatever
they are (it lifts a weak reflection's model amplitude off zero). The comparison is paired, reflection by reflection, $|F_o - |F_\mathrm{null}|| - |F_o - |F_\mathrm{twin}||$
over the free set, its mean against its standard error (`REFINE_TWIN_NULL_Z`): two refinements of one model
differ in R-free by a few thousandths from their trajectories alone, more than a fixed tie, and a fraction of one
or two percent is not resolved by them. A pseudo-symmetry the model already carries (non-crystallographic
**The decision.** The twinned refinement is kept where its fraction is above 0 and the law beats its own **null**
on the free reflections by at least three standard errors - the null being the same twinned model at the same
fraction with every reflection's mate replaced by an unrelated reflection of the same resolution and of the same
set (a free reflection's by a free one, as a law's are): what the averaging of two intensities buys whatever they
are (it lifts a weak reflection's model amplitude off zero). The comparison is paired, reflection by reflection,
$|F_o - |F_\mathrm{null}|| - |F_o - |F_\mathrm{twin}||$ over the free set, summed over each free twin pair (the two
mates share their two $F_\mathrm{calc}$, so they are one term, not two independent ones), its mean against its
standard error (`REFINE_TWIN_NULL_Z`). The R-free of the twinned refinement against the untwinned one's does not
enter: two refinements of one model differ in R-free by a few thousandths from their trajectories alone, and any
twin target lowers R by its averaging, so that comparison favours the twinned one with or without a twin. A pseudo-symmetry the model already carries (non-crystallographic
symmetry parallel to the law) predicts $I(T\mathbf{h}) \approx I(\mathbf{h})$, so the law changes nothing and is not
taken for a twin. A law is tested only where the free set keeps its twin pairs together, and is otherwise reported
taken for a twin.
**The data's own evidence.** A lead of two to three standard errors is enough where the merge's intensity statistics
read twinned as well - the L-test itself below 0.44 (`TWINNING_VERDICT` `INDICATED` with
`ESTIMATED_TWIN_FRACTION_FROM` `L_TEST`; not where only the second moment reads low, nor where a pseudo-translation
leaves the L-test unreadable): two independent pieces of evidence, the data's and the model's, then agree.
`REFINE_TWIN_RULE` says which rule decided. The data-side test is not made a gate - it misses twins a model sees, a
pseudo-merohedral law or a near-perfect twin among them - only further evidence where the null alone is short.
`--mode refine` has no merge statistics, so there the null alone decides.
A law is tested only where the free set keeps its twin pairs together, and is otherwise reported
(`REFINE_TWIN` `FREE_SET_SPLIT`) rather than tested. The free set rugnux draws (§10.7) keeps together the pairs of
every law that is a rotation of the cell's own basis vectors; a law of a centred lattice that is not (its matrix
fractional in the cell's basis) can be split by it, and so can a free set `--mode refine` reads from a file. Of several laws, the one the placed model asks the largest fraction of is refined. Each of the
+3 -2
View File
@@ -1251,8 +1251,9 @@ Where a model was placed and refined ([Refinement](RUGNUX_MODEL_AND_REFERENCE.md
| `REFINE_TWIN_LAW`, `REFINE_TWIN_FRACTION` | The twin law tested (an hkl operator in the data's basis) and the fraction it refined to |
| `REFINE_TWIN_SPLIT_PAIRS` | Where `FREE_SET_SPLIT`: the twin pairs with one member in the working and one in the free set |
| `REFINE_TWIN_R_WORK`, `REFINE_TWIN_R_FREE`, `REFINE_UNTWINNED_R_WORK`, `REFINE_UNTWINNED_R_FREE` | The twinned and the untwinned refinement, each where it ended (its solvent built) |
| `REFINE_TWIN_R_FREE_NULL` | The twinned model's R-free at its fraction with every mate replaced by an unrelated reflection of the same resolution: what averaging alone buys. `REFINE_TWIN_FRACTION_NULL` (*dev*) is the fraction those unrelated mates ask for, about 0 |
| `REFINE_TWIN_NULL_Z` | The law's lead over that null on the free reflections, paired reflection by reflection, in standard errors; `TWINNED` needs 3 |
| `REFINE_TWIN_R_FREE_NULL` | The twinned model's R-free at its fraction with every mate replaced by an unrelated reflection of the same resolution and set: what averaging alone buys, and the like-for-like comparison for `REFINE_TWIN_R_FREE`. `REFINE_TWIN_FRACTION_NULL` (*dev*) is the fraction those unrelated mates ask for, about 0 |
| `REFINE_TWIN_NULL_Z` | The law's lead over that null on the free reflections, paired reflection by reflection (a free twin pair as one term), in standard errors; `TWINNED` needs 3, or 2 where the data read twinned (`REFINE_TWIN_RULE`) |
| `REFINE_TWIN_RULE` | Where `UNTWINNED` or `TWINNED`, the rule that decided ([§15.12](CPU_DATA_ANALYSIS_REFINEMENT.md#1512-twinned-data)): `NULL` (the null alone, at 3 standard errors) or `NULL+LTEST` (taken at 2-3 standard errors because the merge's L-test reads twinned as well: `TWINNING_VERDICT` `INDICATED` with `ESTIMATED_TWIN_FRACTION_FROM` `L_TEST`). In `--mode refine` there are no merge statistics, and it is always `NULL` |
| `REFINE_BLOBS` | The blobs of the refined mFo-DFc map (above 3 σ, peaking above 4 σ) |
| `REFINE_BLOBS_LIGAND`, `REFINE_BLOBS_WATER`, `REFINE_BLOBS_ION`, `REFINE_BLOBS_MODEL_ERROR` | How many of them are read as each ([§15.9](CPU_DATA_ANALYSIS_REFINEMENT.md#159-unexplained-density-blobs)); the rest are unassigned |
| `REFINE_SOLVENT_SITES` | The waters and single-atom ions of the refined model the solvent audit scored ([§15.10](CPU_DATA_ANALYSIS_REFINEMENT.md#1510-the-solvent-audit)) |
@@ -1266,24 +1266,38 @@ void MakeMaps(Context &ctx, const std::vector<MergedReflection> &merged, double
}
// The null of a twinned refinement, for the model as the last Refit() left it: its R-free at the same fraction with
// every reflection's mate replaced by an unrelated observed one of the same resolution - the reflections in order of
// 1/d^2 paired TWIN_NULL_SHIFT apart, beyond any run of equal 1/d^2 a twin law relates. What the averaging of two
// intensities buys on its own, whatever they are.
constexpr size_t TWIN_NULL_SHIFT = 8;
// The twinned refinement is taken only where the law beats its null on the free set by this many standard errors.
// every reflection's mate replaced by an unrelated observed one of the same resolution and the same set (free with
// free, work with work, as a law's twin pairs are). The reflections of each set are put in order of 1/d^2 with the two
// members of a twin pair next to each other, and paired TWIN_NULL_SHIFT apart - past the reflection's own pair, and no
// further, as the free set is sparse and a wider step takes the mate to another resolution, which a null pays for.
// What the averaging of two intensities buys on its own, whatever they are.
constexpr size_t TWIN_NULL_SHIFT = 2;
// The twinned refinement is taken only where the law beats its null on the free set by this many standard errors -
// or by TWIN_NULL_Z_DATA where the data's own intensity statistics read twinned, a second, independent piece of
// evidence (AdoptTwinLaw).
constexpr double TWIN_NULL_Z = 3.0;
constexpr double TWIN_NULL_Z_DATA = 2.0;
void TwinNull(const Context &ctx, RefinementResult &result) {
RefinementReflections null = ctx.refl;
const size_t n = null.Size();
null.mate_hkl.clear();
std::vector<int> order(n);
std::iota(order.begin(), order.end(), 0);
std::stable_sort(order.begin(), order.end(), [&](int a, int b) { return null.inv_d2[a] < null.inv_d2[b]; });
null.mate.assign(n, 0);
for (size_t r = 0; r < n; r++) {
// in blocks of twice the shift, the first half with the second
const size_t p = r % (2 * TWIN_NULL_SHIFT) < TWIN_NULL_SHIFT ? r + TWIN_NULL_SHIFT : r - TWIN_NULL_SHIFT;
null.mate[order[r]] = order[p < n ? p : r];
for (const bool free : {false, true}) {
std::vector<int> order;
for (size_t m = 0; m < n; m++)
if ((null.free[m] != 0) == free)
order.push_back(static_cast<int>(m));
// a twin pair is ordered at the 1/d^2 of its first member
auto pair_of = [&](int m) { return ctx.refl.mate[m] < static_cast<int>(n) ? std::min(m, ctx.refl.mate[m]) : m; };
std::stable_sort(order.begin(), order.end(), [&](int a, int b) {
const int pa = pair_of(a), pb = pair_of(b);
return std::make_pair(null.inv_d2[pa], pa) < std::make_pair(null.inv_d2[pb], pb);
});
for (size_t r = 0; r < order.size(); r++) {
// in blocks of twice the shift, the first half with the second
const size_t p = r % (2 * TWIN_NULL_SHIFT) < TWIN_NULL_SHIFT ? r + TWIN_NULL_SHIFT : r - TWIN_NULL_SHIFT;
null.mate[order[r]] = order[p < order.size() ? p : r];
}
}
// at the law's fraction: averaging any two intensities lifts the weak reflections' model amplitudes off zero
// by itself (the fraction the unrelated mates ask for, reported beside it, is about 0)
@@ -1292,20 +1306,27 @@ void TwinNull(const Context &ctx, RefinementResult &result) {
const std::vector<std::complex<float>> f_null = TwinnedModel(null, rows, ctx.fcalc);
result.twin_r_free_null = ComputeRFactors(null, f_null).r_free;
// How far the law beats its null on the free reflections, paired: the mean of |Fo - |F_null|| - |Fo - |F_twin||
// over its standard error.
// over its standard error. The two members of a free twin pair share their two F_calc, so their terms are summed
// into one: the standard error counts twin pairs, not reflections.
const std::vector<std::complex<float>> f_law = TwinnedModel(ctx.refl, ctx.rows, ctx.fcalc);
double sum = 0, sum2 = 0;
int n_free = 0;
std::vector<double> d(n, 0.0);
for (size_t m = 0; m < n; m++)
if (ctx.refl.free[m]) {
const double d = std::fabs(ctx.refl.fo[m] - std::abs(f_null[m])) - std::fabs(ctx.refl.fo[m] - std::abs(f_law[m]));
sum += d;
sum2 += d * d;
n_free++;
}
if (n_free > 1) {
const double mean = sum / n_free, var = (sum2 - n_free * mean * mean) / (n_free - 1);
result.twin_null_z = var > 0 ? mean / std::sqrt(var / n_free) : 0.0;
if (ctx.refl.free[m])
d[m] = std::fabs(ctx.refl.fo[m] - std::abs(f_null[m])) - std::fabs(ctx.refl.fo[m] - std::abs(f_law[m]));
double sum = 0, sum2 = 0;
int n_pairs = 0;
for (size_t m = 0; m < n; m++) {
const size_t t = ctx.refl.mate[m];
if (!ctx.refl.free[m] || (t < n && t < m))
continue; // the pair was counted at its first member
const double pair = t < n && t != m ? d[m] + d[t] : d[m];
sum += pair;
sum2 += pair * pair;
n_pairs++;
}
if (n_pairs > 1) {
const double mean = sum / n_pairs, var = (sum2 - n_pairs * mean * mean) / (n_pairs - 1);
result.twin_null_z = var > 0 ? mean / std::sqrt(var / n_pairs) : 0.0;
}
}
@@ -3028,11 +3049,12 @@ double FinalRWork(const RefinementResult &r) {
// The twin hypothesis beside the untwinned refinement `result`: the placed model refined again, solvent build and
// displacement model included, against the data as twinned by a law of the lattice, and the twinned refinement taken
// (into `result`) where its R-free is below the untwinned one's and the law beats its own null (TwinNull) on the free
// reflections by TWIN_NULL_Z standard errors - so that neither the drop in R any twin target buys by averaging, nor a
// pseudo-symmetry the model already carries (it predicts I(Th) ~ I(h), and the law then changes nothing) is read as
// a twin. A law is tested only where the free set keeps its twin pairs together; of several, the one the placed
// model asks the largest fraction of.
// (into `result`) where the law beats its own null (TwinNull) on the free reflections by TWIN_NULL_Z standard errors,
// or by TWIN_NULL_Z_DATA where the caller says the data read twinned (AdoptTwinLaw) - not on its R-free against the
// untwinned one's, which compares two refinement paths and so carries their noise - so that neither the drop in R
// any twin target buys by averaging, nor a pseudo-symmetry the model already carries (it predicts I(Th) ~ I(h), and
// the law then changes nothing) is read as a twin. A law is tested only where the free set keeps its twin pairs
// together; of several, the one the placed model asks the largest fraction of.
void DecideTwin(const gemmi::Structure &placed, const std::vector<MergedReflection> &merged,
const gemmi::UnitCell &cell, const gemmi::SpaceGroup &sg, double wavelength_A, Logger &logger,
const RefinementSettings &settings, RefinementResult &result) {
@@ -3091,17 +3113,20 @@ void DecideTwin(const gemmi::Structure &placed, const std::vector<MergedReflecti
result.twin_r_free_null = twinned.twin_r_free_null;
result.twin_fraction_null = twinned.twin_fraction_null;
result.twin_null_z = twinned.twin_null_z;
const bool adopt = twinned.twin_fraction > 0 && result.twin_r_free < result.untwinned_r_free - R_FREE_TIE &&
result.twin_null_z >= TWIN_NULL_Z;
const bool adopt = AdoptTwinLaw(twinned.twin_fraction, result.twin_null_z, settings.twin_data_indicated);
result.twin_decision = twinned.twin_fraction == 0 ? "ZERO_FRACTION" : adopt ? "TWINNED" : "UNTWINNED";
if (twinned.twin_fraction > 0)
result.twin_rule = adopt && result.twin_null_z < TWIN_NULL_Z ? "NULL+LTEST" : "NULL";
logger.Info("Refinement, twin law {}: fraction {:.3f}, R-work {:.4f}, R-free {:.4f}; untwinned R-free {:.4f}; null "
"(mates unrelated, fraction {:.3f}) R-free {:.4f}, the law ahead of it by {:.1f} standard errors -> {}",
"(mates unrelated, fraction {:.3f}) R-free {:.4f}, the law ahead of it by {:.1f} standard errors{} -> {}",
result.twin_law, result.twin_fraction, result.twin_r_work, result.twin_r_free, result.untwinned_r_free,
result.twin_fraction_null, result.twin_r_free_null, result.twin_null_z, result.twin_decision);
result.twin_fraction_null, result.twin_r_free_null, result.twin_null_z,
settings.twin_data_indicated ? ", the data read twinned" : "", result.twin_decision);
if (!adopt)
return;
// the twinned refinement is the result, with the evidence that chose it
twinned.twin_decision = result.twin_decision;
twinned.twin_rule = result.twin_rule;
twinned.twin_laws = result.twin_laws;
twinned.twin_r_work = result.twin_r_work;
twinned.twin_r_free = result.twin_r_free;
@@ -3112,6 +3137,10 @@ void DecideTwin(const gemmi::Structure &placed, const std::vector<MergedReflecti
} // namespace
bool AdoptTwinLaw(double twin_fraction, double null_z, bool data_indicated) {
return twin_fraction > 0 && (null_z >= TWIN_NULL_Z || (data_indicated && null_z >= TWIN_NULL_Z_DATA));
}
RefinementResult RefineModel(const gemmi::Structure &placed, const std::vector<MergedReflection> &merged,
const gemmi::UnitCell &cell, const gemmi::SpaceGroup &sg, double wavelength_A,
Logger &logger, const RefinementSettings &settings) {
@@ -75,6 +75,9 @@ struct RefinementSettings {
// A refinement against the data as twinned by each twin law the lattice admits (RefinementTwinLaws), beside the
// untwinned one, and the model of whichever the evidence favours (RefineModel).
bool twin = true;
// The data's own intensity statistics read twinned - the caller's verdict, from the merge (the L-test): the law
// is then taken from 2 standard errors over its null instead of 3 (AdoptTwinLaw).
bool twin_data_indicated = false;
bool gpu = true;
size_t nthreads = 1;
size_t gpu_resident_bytes = 0; // held on the card by what runs beside
@@ -219,8 +222,15 @@ struct RefinementResult {
double untwinned_r_work = NAN, untwinned_r_free = NAN;
double twin_r_free_null = NAN, twin_fraction_null = NAN;
double twin_null_z = NAN; // the law's lead over the null on the free set, in standard errors
// Where tested to UNTWINNED or TWINNED, the rule that decided: NULL (the null alone, 3 standard errors), or
// NULL+LTEST (taken at 2-3 standard errors because the data read twinned, settings.twin_data_indicated).
std::string twin_rule;
};
// Whether the twinned refinement is taken (RefineModel): a fraction above 0, and the law ahead of its null on the
// free set by 3 standard errors - or by 2 where the data's own intensity statistics read twinned as well.
bool AdoptTwinLaw(double twin_fraction, double null_z, bool data_indicated);
// `placed`: the model in the data's cell and space group. With settings.solvent (and restrained, not score_only),
// the weight trajectory is followed by up to SOLVENT_BUILD_ROUNDS rounds of the solvent build, each the trajectory
// again from the model with its solvent changed, kept where it ends at a lower R-free (Refinement.cpp). `merged`: the
+18 -13
View File
@@ -595,6 +595,8 @@ namespace {
Add(s, KeyReal("REFINE_TWIN_R_FREE_NULL", rf.twin_r_free_null, "{:.4f}"));
Add(s, KeyReal("REFINE_TWIN_FRACTION_NULL", rf.twin_fraction_null, "{:.3f}", true));
Add(s, KeyReal("REFINE_TWIN_NULL_Z", rf.twin_null_z, "{:.1f}"));
if (!rf.twin_rule.empty())
Add(s, KeyEnum("REFINE_TWIN_RULE", rf.twin_rule, {"NULL", "NULL+LTEST"}));
Add(s, KeyReal("REFINE_UNTWINNED_R_WORK", rf.untwinned_r_work, "{:.4f}"));
Add(s, KeyReal("REFINE_UNTWINNED_R_FREE", rf.untwinned_r_free, "{:.4f}"));
}
@@ -672,6 +674,20 @@ namespace {
}
}
// One verdict, in the words every other place uses (TwinningVerdictLine): INDICATED is a twin in a class that
// admits one; SYMMETRY_SUSPECT is the same low reading in a holohedral class, where it says an adopted operator
// averages unequal intensities. NOT_READABLE where <|L|> reads outside the range any twin fraction gives: then
// neither twinning nor a false operator is read from it. NO_TWIN_LAW: twin-like, but the lattice admits no twin law
// for the group (see TwinningAnalysisResult::twin_like_without_twin_law) - a question about the symmetry.
std::string TwinningVerdictKey(const ProcessResult &r) {
const auto &tw = r.twinning;
return TwinStatisticsAbnormal(r) ? "NOT_READABLE"
: tw.adopted_operators_suspect ? "SYMMETRY_SUSPECT"
: tw.twinning_suspected ? "INDICATED"
: tw.twin_like_without_twin_law ? "NO_TWIN_LAW"
: "NO_INDICATION";
}
// ---------------------------------------------------------------------------- the emitter
// Populates the document. Nothing here writes to a stream: the order of the FILE is the renderer's
// business, which is what lets the SUMMARY sit above the sections whose warnings it reports.
@@ -2153,20 +2169,9 @@ ReportDocument BuildReportDocument(const std::string &output_prefix,
: "UNAFFECTED",
{"UNAFFECTED", "REPAIRED", "UNREADABLE"}));
Add(s, KeyReal("SECOND_MOMENT_I", tw.second_moment, "{:.4f}"));
// One verdict, in the words every other place uses (TwinningVerdictLine): INDICATED is a twin
// in a class that admits one; SYMMETRY_SUSPECT is the same low reading in a holohedral class,
// where it says an adopted operator averages unequal intensities.
// NOT_READABLE where <|L|> reads outside the range any twin fraction gives: then neither
// twinning nor a false operator is read from it.
const bool abnormal = TwinStatisticsAbnormal(result);
// NO_TWIN_LAW: twin-like, but the lattice admits no twin law for the group (see
// TwinningAnalysisResult::twin_like_without_twin_law) - a question about the symmetry.
Add(s, KeyEnum("TWINNING_VERDICT",
abnormal ? "NOT_READABLE"
: tw.adopted_operators_suspect ? "SYMMETRY_SUSPECT"
: tw.twinning_suspected ? "INDICATED"
: tw.twin_like_without_twin_law ? "NO_TWIN_LAW"
: "NO_INDICATION",
// One verdict (TwinningVerdictKey), the one refinement's twin decision reads too.
Add(s, KeyEnum("TWINNING_VERDICT", TwinningVerdictKey(result),
{"NO_INDICATION", "INDICATED", "SYMMETRY_SUSPECT", "NO_TWIN_LAW", "NOT_READABLE"}));
// The fraction comes from the statistic that carries the verdict, and the key beside it says
// which. None is quoted where it cannot be measured: merged under a false operator, or from a
+4
View File
@@ -59,6 +59,10 @@ inline constexpr const char *RUGNUX_ACKNOWLEDGEMENT =
// The overall CC_anom above which the summary says the data carry anomalous signal.
inline constexpr double CC_ANOM_SIGNAL = 0.15;
// TWINNING_VERDICT: NO_INDICATION, INDICATED, SYMMETRY_SUSPECT, NO_TWIN_LAW or NOT_READABLE, from the merge's
// intensity statistics (result.twinning, and the space-group search's merge where one was promoted from).
std::string TwinningVerdictKey(const ProcessResult &result);
// The emitter. Populates the result object every rendering reads from; it writes no text itself.
ReportDocument BuildReportDocument(const std::string &output_prefix,
const std::string &input_file,
+13 -2
View File
@@ -12,6 +12,7 @@
#include "../image_analysis/SensorAbsorption.h"
#include "DiagnosticOutput.h"
#include "RugnuxFirstPass.h"
#include "ResultReport.h"
#include <algorithm>
#include <atomic>
@@ -4383,7 +4384,10 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) {
laws += (laws.empty() ? "" : " ; ") + t;
// What each law reads on the merge: the disagreement H of the reflections it
// relates (Yeates (1997) Methods Enzymol. 276, 344-358: a twin of fraction a reads
// (1-2a)/2, so a >= 0.5 - H), their correlation and their intensity-weighted R.
// (1-2a)/2), their correlation and their intensity-weighted R. 0.5 - H is the
// fraction only where nothing but twinning makes the pairs agree: the median over
// the stronger pairs and a pseudo-symmetry along the law both lower H, and
// untwinned crystals read 0.5 - H up to about 0.25 (measured over the battery).
// Report-only, so a user sees which law a twin refinement would need.
gemmi::Mat33 m;
for (int i = 0; i < 3; ++i)
@@ -4395,7 +4399,8 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) {
law_correlation.emplace(expanded_merge(), law_opts);
if (const auto s = (*law_correlation)(m); s.has_value())
result.twin_law_evidence.push_back(fmt::format(
"{}: {} pairs, H {:.3f} (twin fraction {:.2f} or more), R {:.3f}, CC {:.3f}",
"{}: {} pairs, H {:.3f} (twin fraction {:.2f} if only twinning makes the pairs "
"agree; untwinned crystals read up to about 0.25), R {:.3f}, CC {:.3f}",
t, s->n_pairs, s->h_stat, std::max(0.0, 0.5 - s->h_stat), s->r_stat, s->cc));
else
result.twin_law_evidence.push_back(fmt::format("{}: too few pairs to read", t));
@@ -5380,6 +5385,12 @@ bool Rugnux::ScaleMergeAndSymmetry(PipelineLocals &p) {
settings.gpu_resident_bytes = rsm && rsm->OnGPU() ? rsm->DeviceBytes() : 0;
// A cut (the broker's deadline or cancel, Ctrl+C) cuts a running refinement short, keeping what it has
settings.cancel = cancelled_.RunFlag();
// The data-side evidence for a twin law, from this merge: the L-test itself reads twinned - the verdict
// the report gives (TWINNING_VERDICT INDICATED), carried by the L-test (ESTIMATED_TWIN_FRACTION_FROM
// L_TEST), not by a second moment alone or an L-test a pseudo-translation left unreadable. With it the
// refinement takes a law at a smaller lead over its null (AdoptTwinLaw).
settings.twin_data_indicated = TwinningVerdictKey(result) == "INDICATED"
&& result.twinning.twin_fraction_source == TwinFractionSource::LTest;
const PlacedModelOutput p = *placed_model_out; // a copy: a rival taken below replaces *placed_model_out
// In the group RefinableIn chose: the model's own where it has the rotations of the group written -
// the screws are then the model's, as its fit already took them - so a model whose screws the merge
+14
View File
@@ -1540,6 +1540,18 @@ TEST_CASE("RefinementTwin_laws_and_mates", "[Refinement]") {
CHECK(SetTwinLaw(split, t.cell, t.sg, 2.2, t.law) > 0);
}
TEST_CASE("RefinementTwin_rule_with_data_evidence", "[Refinement]") {
// the null alone: 3 standard errors
CHECK(AdoptTwinLaw(0.1, 3.0, false));
CHECK_FALSE(AdoptTwinLaw(0.1, 2.9, false));
// the data read twinned as well: 2 standard errors
CHECK(AdoptTwinLaw(0.1, 2.0, true));
CHECK(AdoptTwinLaw(0.1, 2.5, true));
CHECK_FALSE(AdoptTwinLaw(0.1, 1.9, true));
// no fraction, no twin, whatever the evidence
CHECK_FALSE(AdoptTwinLaw(0.0, 10.0, true));
}
TEST_CASE("RefinementTwin_gradient", "[Refinement]") {
TwinCase t;
gemmi::Structure truth = random_crystal(t.cell, "P 4", 240, 17);
@@ -1636,7 +1648,9 @@ TEST_CASE("Refinement_twin_decided", "[Refinement]") {
CHECK(r.twin_fraction == Catch::Approx(fraction).margin(0.03));
CHECK(r.r_free < r.untwinned_r_free - 0.02);
CHECK(r.twin_r_free < r.twin_r_free_null - 0.02);
CHECK(r.twin_null_z >= 3.0);
} else {
CHECK((r.twin_decision == "ZERO_FRACTION" || r.twin_null_z < 3.0));
CHECK(r.twin_decision != "TWINNED");
CHECK_FALSE(r.twinned);
// the untwinned refinement, as it is without the twin hypothesis