Remove non-helping stills merge/scaling knobs
Trims three opt-in stills parameters that did not improve data quality on the external-reference (PDB R-free) battery and only added code: - --partiality-uncertainty: the (1-p)/p merge-sigma term was null on KR2/OCP/ NmHR/CTX vs their reference structures (and neutral-to-harmful at higher coefficients); removed the flag, setting and CorrectedSigma term. - --stills-modulation: the detector-plane flat-field surface was net-negative on flooded data; removed the flag, setting and MergeOnTheFly::RefineModulation (the rotation modulation in RotationScaleMerge is unaffected). - --min-indexed-fraction: every value other than the 0.20 default collapsed CC1/2; removed the override flag/setter, keeping the fixed 0.20 acceptance floor. Default behaviour is unchanged (all three were off / at their default). Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
This commit is contained in:
Claude Opus 4.8
@@ -28,11 +28,6 @@ int64_t IndexingSettings::GetViableCellMinSpots() const {
|
||||
return viable_cell_min_spots;
|
||||
}
|
||||
|
||||
IndexingSettings &IndexingSettings::MinIndexedSpotFraction(float input) {
|
||||
min_indexed_spot_fraction = input;
|
||||
return *this;
|
||||
}
|
||||
|
||||
float IndexingSettings::GetMinIndexedSpotFraction() const {
|
||||
return min_indexed_spot_fraction;
|
||||
}
|
||||
|
||||
@@ -43,7 +43,6 @@ public:
|
||||
IndexingSettings();
|
||||
|
||||
IndexingSettings& ViableCellMinSpots(int64_t input);
|
||||
IndexingSettings& MinIndexedSpotFraction(float input);
|
||||
IndexingSettings& Algorithm(IndexingAlgorithmEnum input);
|
||||
IndexingSettings& FFT_MaxUnitCell_A(float input);
|
||||
IndexingSettings& FFT_MinUnitCell_A(float input);
|
||||
|
||||
@@ -151,15 +151,6 @@ bool ScalingSettings::GetCorrectionSurfaces() const {
|
||||
return correction_surfaces;
|
||||
}
|
||||
|
||||
ScalingSettings &ScalingSettings::StillsModulation(bool input) {
|
||||
stills_modulation = input;
|
||||
return *this;
|
||||
}
|
||||
|
||||
bool ScalingSettings::GetStillsModulation() const {
|
||||
return stills_modulation;
|
||||
}
|
||||
|
||||
ScalingSettings &ScalingSettings::StillsPartialityRefine(bool input) {
|
||||
stills_partiality_refine = input;
|
||||
return *this;
|
||||
@@ -229,15 +220,6 @@ double ScalingSettings::GetCaptureUncertaintyCoeff() const {
|
||||
return capture_uncertainty_coeff;
|
||||
}
|
||||
|
||||
ScalingSettings &ScalingSettings::PartialityUncertaintyCoeff(double input) {
|
||||
partiality_uncertainty_coeff = input;
|
||||
return *this;
|
||||
}
|
||||
|
||||
double ScalingSettings::GetPartialityUncertaintyCoeff() const {
|
||||
return partiality_uncertainty_coeff;
|
||||
}
|
||||
|
||||
ScalingSettings &ScalingSettings::MinCapturedFraction(double input) {
|
||||
if (input < 0.0)
|
||||
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
||||
|
||||
@@ -27,7 +27,6 @@ class ScalingSettings {
|
||||
// a fraction f<1 of its rocking curve is extrapolated, and the unobserved (1-f) carries a
|
||||
// systematic error ~coeff*(1-f)*I that plain counting sigma misses. 0 = off (baseline).
|
||||
double capture_uncertainty_coeff = 0.0;
|
||||
double partiality_uncertainty_coeff = 0.0;
|
||||
// Full-level captured-fraction floor for the rot3d combine: drop a reconstructed full whose rocking
|
||||
// curve was only fractionally captured (sum of its partials' partiality < this). Distinct from
|
||||
// min_partiality, which gates individual partials; this gates the assembled full. 0 = off (baseline).
|
||||
@@ -46,10 +45,6 @@ class ScalingSettings {
|
||||
bool correction_surfaces = true;
|
||||
// Absorption-surface refinement iteration count (used when correction_surfaces is on).
|
||||
int absorption_iter = 3;
|
||||
// Detector-plane modulation (flat-field) correction for the STILLS merge path (ScaleOnTheFly +
|
||||
// MergeOnTheFly), which otherwise has no correction surfaces. Off by default (experimental; the
|
||||
// rotation path fits its own modulation via RotationScaleMerge). Enabled by rugnux --stills-modulation.
|
||||
bool stills_modulation = false;
|
||||
|
||||
// Physical partiality post-refinement for the STILLS merge (StillsPartialityRefine): refine a per-crystal
|
||||
// orientation tilt against the running merge, recompute each reflection's partiality from the refined
|
||||
@@ -100,14 +95,12 @@ public:
|
||||
ScalingSettings& MinPartiality(double min_partiality);
|
||||
ScalingSettings& ForcedMosaicity(std::optional<double> input);
|
||||
ScalingSettings& CaptureUncertaintyCoeff(double input);
|
||||
ScalingSettings& PartialityUncertaintyCoeff(double input);
|
||||
ScalingSettings& MinCapturedFraction(double input);
|
||||
ScalingSettings& MinCCForImage(double min_cc_for_image);
|
||||
ScalingSettings& OutlierRejectNsigma(double input);
|
||||
ScalingSettings& ScaleFulls(bool input);
|
||||
ScalingSettings& AbsorptionIter(int input);
|
||||
ScalingSettings& CorrectionSurfaces(bool input);
|
||||
ScalingSettings& StillsModulation(bool input);
|
||||
ScalingSettings& StillsPartialityRefine(bool input);
|
||||
ScalingSettings& ExpectedVarianceMerge(bool input);
|
||||
ScalingSettings& SmoothGDegrees(double input);
|
||||
@@ -141,14 +134,12 @@ public:
|
||||
[[nodiscard]] double GetMinPartiality() const;
|
||||
[[nodiscard]] std::optional<double> GetForcedMosaicity() const;
|
||||
[[nodiscard]] double GetCaptureUncertaintyCoeff() const;
|
||||
[[nodiscard]] double GetPartialityUncertaintyCoeff() const;
|
||||
[[nodiscard]] double GetMinCapturedFraction() const;
|
||||
[[nodiscard]] double GetMinCCForImage() const;
|
||||
[[nodiscard]] double GetOutlierRejectNsigma() const;
|
||||
[[nodiscard]] bool GetScaleFulls() const;
|
||||
[[nodiscard]] int GetAbsorptionIter() const;
|
||||
[[nodiscard]] bool GetCorrectionSurfaces() const;
|
||||
[[nodiscard]] bool GetStillsModulation() const;
|
||||
[[nodiscard]] bool GetStillsPartialityRefine() const;
|
||||
[[nodiscard]] bool GetExpectedVarianceMerge() const;
|
||||
[[nodiscard]] double GetSmoothGDegrees() const;
|
||||
|
||||
@@ -93,7 +93,7 @@ void MergeOnTheFly::AddImage(const IntegrationOutcome &outcome, int64_t image_id
|
||||
continue;
|
||||
auto hkl = generator(r);
|
||||
auto hkl_key = hkl.pack();
|
||||
sigma_corr = CorrectedSigma(I_corr, sigma_corr, r.image_scale_corr, hkl_key, r.partiality);
|
||||
sigma_corr = CorrectedSigma(I_corr, sigma_corr, r.image_scale_corr, hkl_key);
|
||||
|
||||
// Robust outlier rejection: drop this observation if it sits more than
|
||||
// reject_nsigma error-model sigmas from the reflection's median. Needs the active
|
||||
@@ -131,139 +131,8 @@ void MergeOnTheFly::AddImage(const IntegrationOutcome &outcome, int64_t image_id
|
||||
}
|
||||
}
|
||||
|
||||
double MergeOnTheFly::RefineModulation(std::vector<IntegrationOutcome> &outcomes) {
|
||||
// A minimum held-out generalizing gain (fraction of the held-out scatter) before the surface is
|
||||
// applied - a margin, so a noise-level "improvement" never engages the correction. Matches the
|
||||
// rotation ApplyCellSurface gate.
|
||||
constexpr double CV_MIN_RELATIVE_GAIN = 0.02;
|
||||
constexpr int NB = 16;
|
||||
const int ncell = NB * NB;
|
||||
|
||||
// One scaled observation reduced to what the surface fit needs. `parity` (image index & 1) drives the
|
||||
// even/odd cross-validation split - the stills analogue of the rotation frame parity.
|
||||
struct MObs { double I, sigma, corr; float px, py; int32_t group; int parity; int cell; };
|
||||
std::vector<MObs> obs;
|
||||
|
||||
// Accept exactly what AddImage merges (systematic absence, scale/resolution/ice/partiality filters),
|
||||
// and additionally require a finite detector position. Returns the dense ASU-group id or -1.
|
||||
std::unordered_map<uint64_t, int> group_of;
|
||||
auto accept = [&](const Reflection &r, MObs &out) -> bool {
|
||||
if (generator.IsSystematicallyAbsent(r)) return false;
|
||||
if (r.image_scale_corr <= 0.0f || !std::isfinite(r.image_scale_corr)) return false;
|
||||
if (!AcceptReflection(r, high_resolution_limit)) return false;
|
||||
if (exclude_ice_rings && r.on_ice_ring) return false;
|
||||
if (IsMaskedRing(r)) return false;
|
||||
if (r.partiality < min_partiality) return false;
|
||||
if (!std::isfinite(r.predicted_x) || !std::isfinite(r.predicted_y)) return false;
|
||||
const float I_corr = r.I * r.image_scale_corr, sigma_corr = r.sigma * r.image_scale_corr;
|
||||
if (!std::isfinite(I_corr) || !std::isfinite(sigma_corr) || sigma_corr <= 0.0f) return false;
|
||||
const uint64_t key = generator(r).pack();
|
||||
auto [it, inserted] = group_of.try_emplace(key, static_cast<int>(group_of.size()));
|
||||
out.I = r.I; out.sigma = r.sigma; out.corr = r.image_scale_corr;
|
||||
out.px = r.predicted_x; out.py = r.predicted_y; out.group = it->second;
|
||||
return true;
|
||||
};
|
||||
|
||||
float pxmin = std::numeric_limits<float>::infinity(), pxmax = -pxmin, pymin = pxmin, pymax = -pxmin;
|
||||
for (size_t i = 0; i < outcomes.size(); ++i) {
|
||||
const int parity = static_cast<int>(i & 1);
|
||||
for (const auto &r : outcomes[i].reflections) {
|
||||
MObs m{};
|
||||
if (!accept(r, m)) continue;
|
||||
m.parity = parity;
|
||||
pxmin = std::min(pxmin, m.px); pxmax = std::max(pxmax, m.px);
|
||||
pymin = std::min(pymin, m.py); pymax = std::max(pymax, m.py);
|
||||
obs.push_back(m);
|
||||
}
|
||||
}
|
||||
const int n_groups = static_cast<int>(group_of.size());
|
||||
if (obs.size() < static_cast<size_t>(8 * ncell) || !(pxmax > pxmin) || !(pymax > pymin))
|
||||
return 0.0; // too sparse to over-determine a 16x16 surface, or degenerate detector footprint
|
||||
|
||||
const float sx = NB / (pxmax - pxmin), sy = NB / (pymax - pymin);
|
||||
for (auto &m : obs) {
|
||||
const int ix = std::clamp(static_cast<int>((m.px - pxmin) * sx), 0, NB - 1);
|
||||
const int iy = std::clamp(static_cast<int>((m.py - pymin) * sy), 0, NB - 1);
|
||||
m.cell = ix * NB + iy;
|
||||
}
|
||||
|
||||
// Fit the per-cell factor over {parity subset} (parity < 0 = all obs), n_iter alternating rounds against
|
||||
// that subset's own inverse-variance reference: Tikhonov pull to 1, gauge-fixed to a den-weighted
|
||||
// geometric mean of 1 so it never drifts the overall scale. Mirrors RotationScaleMerge::ApplyCellSurface.
|
||||
constexpr int N_ITER = 3;
|
||||
auto fit_surface = [&](int parity) -> std::vector<double> {
|
||||
std::vector<double> A(ncell, 1.0);
|
||||
for (int it = 0; it < N_ITER; ++it) {
|
||||
std::vector<double> sw(n_groups, 0.0), swI(n_groups, 0.0);
|
||||
for (const auto &o : obs) {
|
||||
if (parity >= 0 && o.parity != parity) continue;
|
||||
const double a = A[o.cell], sc = o.sigma * o.corr * a, w = 1.0 / (sc * sc);
|
||||
sw[o.group] += w; swI[o.group] += w * o.I * o.corr * a;
|
||||
}
|
||||
std::vector<double> num(ncell, 0.0), den(ncell, 0.0);
|
||||
for (const auto &o : obs) {
|
||||
if ((parity >= 0 && o.parity != parity) || sw[o.group] <= 0.0) continue;
|
||||
const double Iref = swI[o.group] / sw[o.group], a = A[o.cell];
|
||||
const double Is = o.I * o.corr * a, sc = o.sigma * o.corr * a;
|
||||
if (!std::isfinite(Iref) || Iref <= 0.0 || !(Is > 0.0) || !(sc > 0.0)) continue;
|
||||
const double w = 1.0 / (sc * sc);
|
||||
num[o.cell] += w * Is * Iref; den[o.cell] += w * Is * Is;
|
||||
}
|
||||
std::vector<double> dsorted = den;
|
||||
std::nth_element(dsorted.begin(), dsorted.begin() + dsorted.size() / 2, dsorted.end());
|
||||
const double lambda = 0.1 * std::max(1e-30, dsorted[dsorted.size() / 2]);
|
||||
double logsum = 0.0, wsum = 0.0;
|
||||
std::vector<double> upd(ncell, 1.0);
|
||||
for (int c = 0; c < ncell; ++c) upd[c] = (num[c] + lambda) / (den[c] + lambda);
|
||||
for (int c = 0; c < ncell; ++c) if (den[c] > 0.0) { logsum += den[c] * std::log(upd[c]); wsum += den[c]; }
|
||||
const double gm = wsum > 0.0 ? std::exp(logsum / wsum) : 1.0;
|
||||
for (int c = 0; c < ncell; ++c) A[c] = std::clamp(A[c] * upd[c] / gm, 0.25, 4.0);
|
||||
}
|
||||
return A;
|
||||
};
|
||||
// Sigma-independent (R-meas-like) agreement of the held-out equivalents: sum|Is - Iref| / sum|Iref|.
|
||||
// A fractional metric cannot be gamed by a surface that merely reshapes sigma via corr.
|
||||
auto score = [&](int parity, const std::vector<double> &A) -> double {
|
||||
std::vector<double> sw(n_groups, 0.0), swI(n_groups, 0.0);
|
||||
for (const auto &o : obs) {
|
||||
if (o.parity != parity) continue;
|
||||
const double a = A[o.cell], Is = o.I * o.corr * a, sc = o.sigma * o.corr * a, w = 1.0 / (sc * sc);
|
||||
sw[o.group] += w; swI[o.group] += w * Is;
|
||||
}
|
||||
double num = 0.0, den = 0.0;
|
||||
for (const auto &o : obs) {
|
||||
if (o.parity != parity || sw[o.group] <= 0.0) continue;
|
||||
const double a = A[o.cell], Is = o.I * o.corr * a, Iref = swI[o.group] / sw[o.group];
|
||||
if (!std::isfinite(Iref) || Iref <= 0.0) continue;
|
||||
num += std::abs(Is - Iref); den += Iref;
|
||||
}
|
||||
return den > 0.0 ? num / den : 0.0;
|
||||
};
|
||||
|
||||
// Cross-validate: fit on even images, score the held-out odd equivalents (and vice versa). Apply the
|
||||
// full-data surface only if the held-out agreement improves by a clear margin.
|
||||
const std::vector<double> ident(ncell, 1.0);
|
||||
const std::vector<double> A_even = fit_surface(0), A_odd = fit_surface(1);
|
||||
const double base = score(1, ident) + score(0, ident);
|
||||
const double gain = base - (score(1, A_even) + score(0, A_odd));
|
||||
if (!(gain > CV_MIN_RELATIVE_GAIN * base))
|
||||
return 0.0; // not cross-validated: the correction stays a no-op (caller logs)
|
||||
const std::vector<double> A = fit_surface(-1);
|
||||
|
||||
// Fold the surface into each accepted reflection's image_scale_corr (recompute its cell deterministically).
|
||||
for (auto &outcome : outcomes)
|
||||
for (auto &r : outcome.reflections) {
|
||||
MObs m{};
|
||||
if (!accept(r, m)) continue;
|
||||
const int ix = std::clamp(static_cast<int>((m.px - pxmin) * sx), 0, NB - 1);
|
||||
const int iy = std::clamp(static_cast<int>((m.py - pymin) * sy), 0, NB - 1);
|
||||
r.image_scale_corr = static_cast<float>(r.image_scale_corr * A[ix * NB + iy]);
|
||||
}
|
||||
return gain / std::max(base, 1e-30); // held-out gain fraction (caller logs)
|
||||
}
|
||||
|
||||
float MergeOnTheFly::CorrectedSigma(float I_corr, float sigma_corr, float image_scale_corr,
|
||||
uint64_t hkl_key, float partiality) const {
|
||||
uint64_t hkl_key) const {
|
||||
if (!error_model_active)
|
||||
return sigma_corr;
|
||||
|
||||
@@ -288,21 +157,8 @@ float MergeOnTheFly::CorrectedSigma(float I_corr, float sigma_corr, float image_
|
||||
a_var = base;
|
||||
}
|
||||
|
||||
double v = error_model_a * a_var
|
||||
+ (error_model_b * I_for_b) * (error_model_b * I_for_b);
|
||||
|
||||
// Partiality-model uncertainty: a reflection recorded at fraction p carries a systematic intensity
|
||||
// error ~ (dp/p) that is proportional to <I> and grows as p falls - plain counting sigma misses it,
|
||||
// so strong low-p partials would otherwise be over-trusted. This is the stills-partiality analog of
|
||||
// the rotation --capture-uncertainty term ((1-captured_fraction)*I in RotationScaleMerge). Inert when
|
||||
// partiality == 1 (no stills partiality model). Gated on a real systematic (error_model_b > 1, i.e.
|
||||
// ISa < 1): on weak counting-limited data (small b) it would only over-concentrate the merge and hurt.
|
||||
const double c = scaling_settings.GetPartialityUncertaintyCoeff();
|
||||
if (c > 0.0 && error_model_b > 1.0) {
|
||||
const double one_minus_p = std::max(0.0, std::min(1.0, 1.0 - static_cast<double>(partiality)));
|
||||
const double t = c * I_for_b * one_minus_p;
|
||||
v += t * t;
|
||||
}
|
||||
const double v = error_model_a * a_var
|
||||
+ (error_model_b * I_for_b) * (error_model_b * I_for_b);
|
||||
return (v > 0.0) ? static_cast<float>(std::sqrt(v)) : sigma_corr;
|
||||
}
|
||||
|
||||
|
||||
@@ -130,7 +130,7 @@ class MergeOnTheFly {
|
||||
// using the per-observation I_i instead would over-weight down-fluctuated points.
|
||||
std::unordered_map<uint64_t, float> error_model_mean_I;
|
||||
[[nodiscard]] float CorrectedSigma(float I_corr, float sigma_corr, float image_scale_corr,
|
||||
uint64_t hkl_key, float partiality) const;
|
||||
uint64_t hkl_key) const;
|
||||
|
||||
// Optional per-observation outlier rejection: drop observations whose corrected
|
||||
// intensity lies more than reject_nsigma error-model sigmas from the reflection's
|
||||
@@ -168,16 +168,6 @@ public:
|
||||
// order (or threading) of AddImage calls - not a draw from a shared RNG in call order.
|
||||
void AddImage(const IntegrationOutcome& outcome, int64_t image_id, bool cc_mask = false);
|
||||
|
||||
// Detector-plane modulation (flat-field) correction for the STILLS path: fit a smooth multiplicative
|
||||
// factor over where a reflection lands on the detector (predicted x,y) against the merged reference and
|
||||
// fold it into each reflection's image_scale_corr, so the following AddImage merge (and the error model)
|
||||
// see the corrected scale. A 16x16 grid, cross-validated (fit even images, score the held-out odd
|
||||
// equivalents by a sigma-independent R-meas-like metric) so it is a no-op when the systematic is absent
|
||||
// or the data too sparse. The stills analogue of RotationScaleMerge::RefineModulation - the stills
|
||||
// ScaleOnTheFly/MergeOnTheFly path otherwise has no correction surfaces. Mutates `outcomes`; call after
|
||||
// scaling and before RefineErrorModel/AddImage. Returns the applied held-out gain fraction (0 = no-op).
|
||||
double RefineModulation(std::vector<IntegrationOutcome> &outcomes);
|
||||
|
||||
// d_min_override, when set, is the effective high-resolution limit for the shell table (used for
|
||||
// the automatic resolution cutoff computed by the caller); otherwise the manual
|
||||
// ScalingSettings high-resolution limit stands. The number of shells is ScalingSettings::ReportShellCount.
|
||||
|
||||
@@ -991,17 +991,6 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b
|
||||
merge_engine.MaskIceRings(masked_ice_rings, config_.spot_finding.ice_ring_width_Q_recipA);
|
||||
if (result.consensus_cell.has_value())
|
||||
merge_engine.ReferenceCell(*result.consensus_cell);
|
||||
// Optional detector-plane modulation (flat-field) correction for stills, folded into each
|
||||
// reflection's scale before the error model and merge (never for the P1 search pass). The
|
||||
// rotation path fits its own modulation inside RotationScaleMerge.
|
||||
if (experiment_.GetScalingSettings().GetStillsModulation() && !for_search) {
|
||||
const double mod_gain = merge_engine.RefineModulation(indexer->GetIntegrationOutcome());
|
||||
if (mod_gain > 0.0)
|
||||
logger.Info("Stills modulation: detector-frame 16x16 surface applied "
|
||||
"(cross-validated, held-out gain {:.1f}%)", 100.0 * mod_gain);
|
||||
else
|
||||
logger.Info("Stills modulation: no cross-validated gain (skipped)");
|
||||
}
|
||||
merge_engine.RefineErrorModel(merge_input);
|
||||
if (merge_engine.ErrorModelActive())
|
||||
logger.Info("Error model: a={:.3f} b={:.3f} ISa={:.1f} chi2={:.2f}", merge_engine.ErrorModelA(),
|
||||
|
||||
@@ -149,10 +149,6 @@ std::string RugnuxCommandLine(const ProcessConfig &config,
|
||||
args.emplace_back("-A");
|
||||
if (sc.GetRefineB())
|
||||
args.emplace_back("-B");
|
||||
if (sc.GetPartialityUncertaintyCoeff() > 0.0)
|
||||
add("--partiality-uncertainty", num(sc.GetPartialityUncertaintyCoeff()));
|
||||
if (sc.GetStillsModulation())
|
||||
args.emplace_back("--stills-modulation");
|
||||
if (!sc.GetStillsPartialityRefine())
|
||||
args.emplace_back("--simple-stills");
|
||||
if (!sc.GetExpectedVarianceMerge())
|
||||
|
||||
@@ -67,7 +67,6 @@ void print_usage() {
|
||||
std::cout << " --min-pix-per-spot <num> Minimum connected strong pixels per spot (default: 2; serial data can index better with 1 + a higher --spot-threshold)" << std::endl;
|
||||
std::cout << " --adaptive-spots Self-calibrating detection: replace the fixed --spot-threshold with a per-resolution-ring threshold set from each image's own noise, so one setting adapts across datasets (no per-dataset --spot-threshold/--spot-sigma tuning)" << std::endl;
|
||||
std::cout << " --spot-false-pixels <num> Adaptive detection operating point: expected noise pixels tolerated per frame (default: 100; implies --adaptive-spots)" << std::endl;
|
||||
std::cout << " --min-indexed-fraction <f> Minimum fraction of in-resolution spots a lattice must index to be accepted (default 0.20); lower to admit weaker/sparser crystals (pair with --min-image-cc to gate the extra spurious ones)" << std::endl;
|
||||
std::cout << " --spot-high-resolution <num> High resolution limit for spot finding (default: 1.5)" << std::endl;
|
||||
std::cout << " --spot-low-resolution <num> Low resolution limit for spot finding, in A (default: 50; lower it, e.g. 24, to exclude the direct-beam halo on weakly-diffracting serial data)" << std::endl;
|
||||
std::cout << " --max-spots <num> Max spot count (default: 250)" << std::endl;
|
||||
@@ -97,7 +96,6 @@ void print_usage() {
|
||||
std::cout << " --smooth-g[=deg] rot3d: smooth per-frame scale G over a deg-degree rotation range (XDS DELPHI-like) before the combine (default: 5 for rot3d; 0 = off)" << std::endl;
|
||||
std::cout << " --relative-b[=deg] rot3d: fit a per-batch relative-B (beyond the single decay slope) over deg-degree batches; cross-validated (default: 10 deg when bare; off otherwise)" << std::endl;
|
||||
std::cout << " --no-scaling-corrections rot3d: disable the (default-on) decay + absorption correction surfaces fitted on the fulls after scale-fulls" << std::endl;
|
||||
std::cout << " --stills-modulation stills: fit a detector-plane modulation (flat-field) surface over the merged reflections (cross-validated; experimental, default off)" << std::endl;
|
||||
std::cout << " --no-expected-variance-merge stills: disable the default expected-variance merge weighting (which rebuilds each weak observation's signal variance at the reflection mean to de-bias the inverse-variance merge); restores observed-sigma weighting" << std::endl;
|
||||
std::cout << " -A, --anomalous Anomalous mode (don't merge Friedel pairs)" << std::endl;
|
||||
std::cout << " -B, --refine-bfactor Refine per image B-factor (stills only)" << std::endl;
|
||||
@@ -107,7 +105,6 @@ void print_usage() {
|
||||
std::cout << " --resolution-shells <num> Number of resolution shells in the reported statistics table (default: 10)" << std::endl;
|
||||
std::cout << " --min-partiality <num> Minimum partiality to accept reflection (default: 0.02)" << std::endl;
|
||||
std::cout << " --capture-uncertainty <num> rot3d: systematic sigma ~num*(1-captured_fraction)*I on under-captured fulls (default: 1.0 for rot3d, 0 otherwise)" << std::endl;
|
||||
std::cout << " --partiality-uncertainty <num> stills: extra merge sigma ~num*(1-partiality)*<I> on partials (auto-gated to error-model b>1 / ISa<1; default 0, ~2.5 recommended)" << std::endl;
|
||||
std::cout << " --min-captured-fraction <num> rot3d: drop a combined full whose rocking curve was captured below this fraction (edge-of-sweep truncated fulls) (default: 0.7 for rotation, 0 otherwise; 0 = off)" << std::endl;
|
||||
std::cout << " --mosaicity <num> Diagnostic: fix the scaling mosaicity (deg) instead of the per-image seed" << std::endl;
|
||||
std::cout << " --reject-outliers <num> Per-observation merge outlier rejection, N sigma from the per-reflection median (default: 6 for rot3d, XDS/DIALS-style; 0 = off)" << std::endl;
|
||||
@@ -148,7 +145,6 @@ enum {
|
||||
OPT_MIN_PIX_PER_SPOT,
|
||||
OPT_ADAPTIVE_SPOTS,
|
||||
OPT_SPOT_FALSE_PIXELS,
|
||||
OPT_MIN_INDEXED_FRACTION,
|
||||
OPT_SPOT_RESOLUTION,
|
||||
OPT_SPOT_LOW_RESOLUTION,
|
||||
OPT_MAX_SPOTS,
|
||||
@@ -175,13 +171,11 @@ enum {
|
||||
OPT_SIMPLE_STILLS,
|
||||
OPT_SCALE_FULLS,
|
||||
OPT_CAPTURE_UNCERTAINTY,
|
||||
OPT_PARTIALITY_UNCERTAINTY,
|
||||
OPT_MIN_CAPTURED_FRACTION,
|
||||
OPT_MOSAICITY,
|
||||
OPT_SMOOTH_G,
|
||||
OPT_RELATIVE_B,
|
||||
OPT_NO_SCALING_CORRECTIONS,
|
||||
OPT_STILLS_MODULATION,
|
||||
OPT_NO_EXPECTED_VARIANCE_MERGE,
|
||||
OPT_DETECT_ICE_RINGS,
|
||||
OPT_NO_SCALE_FULLS,
|
||||
@@ -229,7 +223,6 @@ static option long_options[] = {
|
||||
{"smooth-g", optional_argument, nullptr, OPT_SMOOTH_G},
|
||||
{"relative-b", optional_argument, nullptr, OPT_RELATIVE_B},
|
||||
{"no-scaling-corrections", no_argument, nullptr, OPT_NO_SCALING_CORRECTIONS},
|
||||
{"stills-modulation", no_argument, nullptr, OPT_STILLS_MODULATION},
|
||||
{"no-expected-variance-merge", no_argument, nullptr, OPT_NO_EXPECTED_VARIANCE_MERGE},
|
||||
{"refine", required_argument, nullptr, 'r'},
|
||||
|
||||
@@ -260,13 +253,11 @@ static option long_options[] = {
|
||||
{"min-pix-per-spot", required_argument, nullptr, OPT_MIN_PIX_PER_SPOT},
|
||||
{"adaptive-spots", no_argument, nullptr, OPT_ADAPTIVE_SPOTS},
|
||||
{"spot-false-pixels", required_argument, nullptr, OPT_SPOT_FALSE_PIXELS},
|
||||
{"min-indexed-fraction", required_argument, nullptr, OPT_MIN_INDEXED_FRACTION},
|
||||
{"spot-high-resolution", required_argument, nullptr, OPT_SPOT_RESOLUTION},
|
||||
{"spot-low-resolution", required_argument, nullptr, OPT_SPOT_LOW_RESOLUTION},
|
||||
{"max-spots", required_argument, nullptr, OPT_MAX_SPOTS},
|
||||
{"min-partiality", required_argument, nullptr, OPT_MIN_PARTIALITY},
|
||||
{"capture-uncertainty", required_argument, nullptr, OPT_CAPTURE_UNCERTAINTY},
|
||||
{"partiality-uncertainty", required_argument, nullptr, OPT_PARTIALITY_UNCERTAINTY},
|
||||
{"min-captured-fraction", required_argument, nullptr, OPT_MIN_CAPTURED_FRACTION},
|
||||
{"mosaicity", required_argument, nullptr, OPT_MOSAICITY},
|
||||
{"min-image-cc", required_argument, nullptr, OPT_MIN_IMAGE_CC},
|
||||
@@ -515,7 +506,6 @@ int main(int argc, char **argv) {
|
||||
std::optional<double> smooth_g_deg_arg; // --smooth-g[=deg]; default 5 deg for rot3d, 0 (off) otherwise
|
||||
std::optional<double> relative_b_deg_arg; // --relative-b[=deg]; per-batch relative-B width, 0 (off) unless given
|
||||
bool no_scaling_corrections = false; // --no-scaling-corrections: disable rot3d decay+absorption surfaces
|
||||
bool stills_modulation_flag = false; // --stills-modulation: detector-plane flat-field surface for stills
|
||||
bool no_expected_variance_merge = false; // --no-expected-variance-merge: restore observed-sigma stills merge weighting
|
||||
bool anomalous_mode = false;
|
||||
std::optional<int64_t> space_group_number;
|
||||
@@ -526,7 +516,6 @@ int main(int argc, char **argv) {
|
||||
int64_t min_pix_per_spot = 2;
|
||||
bool adaptive_spots = false;
|
||||
float false_pixels_per_frame = 100.0f;
|
||||
std::optional<float> min_indexed_fraction;
|
||||
bool refine_bfactor = false;
|
||||
std::string ref_mtz;
|
||||
std::string ref_column;
|
||||
@@ -535,7 +524,6 @@ int main(int argc, char **argv) {
|
||||
double min_partiality = 0.02;
|
||||
std::optional<double> min_captured_fraction_arg; // explicit --min-captured-fraction; default depends on rotation
|
||||
std::optional<double> capture_uncertainty_arg; // explicit --capture-uncertainty; default depends on rot3d
|
||||
std::optional<double> partiality_uncertainty_arg; // --partiality-uncertainty (stills partiality merge sigma)
|
||||
std::optional<double> forced_mosaicity_arg; // diagnostic: fix the scaling mosaicity (deg) instead of the per-image seed
|
||||
double min_image_cc = 0.0;
|
||||
int64_t scaling_iter = 3;
|
||||
@@ -765,10 +753,6 @@ int main(int argc, char **argv) {
|
||||
adaptive_spots = true;
|
||||
logger.Info("Adaptive spot detection: expected false pixels/frame set to {:.0f}", false_pixels_per_frame);
|
||||
break;
|
||||
case OPT_MIN_INDEXED_FRACTION:
|
||||
min_indexed_fraction = parse_number_arg<float>(optarg, "--min-indexed-fraction", logger, 0.0f);
|
||||
logger.Info("Minimum indexed-spot fraction for acceptance set to {:.2f}", min_indexed_fraction.value());
|
||||
break;
|
||||
case OPT_SPOT_LOW_RESOLUTION:
|
||||
d_max_spot_finding = parse_number_arg<float>(optarg, "--spot-low-resolution", logger, 0.0f);
|
||||
logger.Info("Low resolution limit for spot finding set to {:.1f} A", d_max_spot_finding);
|
||||
@@ -815,9 +799,6 @@ int main(int argc, char **argv) {
|
||||
case OPT_RELATIVE_B:
|
||||
relative_b_deg_arg = optarg ? parse_double_arg(optarg, "--relative-b", logger) : RELATIVE_B_DEFAULT_DEG;
|
||||
break;
|
||||
case OPT_STILLS_MODULATION:
|
||||
stills_modulation_flag = true;
|
||||
break;
|
||||
case OPT_NO_EXPECTED_VARIANCE_MERGE:
|
||||
no_expected_variance_merge = true;
|
||||
break;
|
||||
@@ -830,9 +811,6 @@ int main(int argc, char **argv) {
|
||||
case OPT_CAPTURE_UNCERTAINTY:
|
||||
capture_uncertainty_arg = parse_double_arg(optarg, "--capture-uncertainty", logger);
|
||||
break;
|
||||
case OPT_PARTIALITY_UNCERTAINTY:
|
||||
partiality_uncertainty_arg = parse_double_arg(optarg, "--partiality-uncertainty", logger);
|
||||
break;
|
||||
case OPT_MIN_CAPTURED_FRACTION:
|
||||
min_captured_fraction_arg = parse_double_arg(optarg, "--min-captured-fraction", logger);
|
||||
break;
|
||||
@@ -1099,7 +1077,6 @@ int main(int argc, char **argv) {
|
||||
scaling_settings.MinCapturedFraction(min_captured_fraction_arg.value_or(
|
||||
(experiment.GetGoniometer().has_value() && !force_still) ? 0.7 : 0.0));
|
||||
scaling_settings.MinCCForImage(min_image_cc / 100.0); // --min-image-cc is percent; the setting is a fraction
|
||||
scaling_settings.StillsModulation(stills_modulation_flag);
|
||||
scaling_settings.StillsPartialityRefine(!simple_stills_flag);
|
||||
scaling_settings.ExpectedVarianceMerge(!no_expected_variance_merge);
|
||||
scaling_settings.OutlierRejectNsigma(
|
||||
@@ -1167,16 +1144,6 @@ int main(int argc, char **argv) {
|
||||
}
|
||||
MergeOnTheFly merge_engine(experiment);
|
||||
merge_engine.ReferenceCell(experiment.GetUnitCell());
|
||||
// Optional detector-plane modulation (flat-field) correction, folded into each reflection's
|
||||
// scale before the error model + merge (mirrors the full-analysis stills path in Rugnux.cpp).
|
||||
if (experiment.GetScalingSettings().GetStillsModulation()) {
|
||||
const double mod_gain = merge_engine.RefineModulation(reflections);
|
||||
if (mod_gain > 0.0)
|
||||
logger.Info("Stills modulation: detector-frame 16x16 surface applied "
|
||||
"(cross-validated, held-out gain {:.1f}%)", 100.0 * mod_gain);
|
||||
else
|
||||
logger.Info("Stills modulation: no cross-validated gain (skipped)");
|
||||
}
|
||||
// Fit the (a, b) error model from symmetry-mate scatter before merging, exactly as the full
|
||||
// pipeline does (Rugnux.cpp). Without this the offline --scale merge would use the identity
|
||||
// model and produce much worse stills intensities (no (b*I)^2 systematic term, no sigma floor).
|
||||
@@ -1404,8 +1371,6 @@ int main(int argc, char **argv) {
|
||||
if (rotation_indexing_range.has_value())
|
||||
indexing_settings.RotationIndexingMinAngularRange_deg(rotation_indexing_range.value());
|
||||
indexing_settings.GeomRefinementAlgorithm(refinement_algorithm);
|
||||
if (min_indexed_fraction.has_value())
|
||||
indexing_settings.MinIndexedSpotFraction(min_indexed_fraction.value());
|
||||
experiment.ImportIndexingSettings(indexing_settings);
|
||||
|
||||
// --detect-ice-rings[=on|off] overrides the value carried in from the dataset (HDF5MetadataSource
|
||||
@@ -1423,7 +1388,6 @@ int main(int argc, char **argv) {
|
||||
scaling_settings.RelativeBDegrees(relative_b_deg_arg.value_or(0.0)); // opt-in only; default off
|
||||
if (no_scaling_corrections)
|
||||
scaling_settings.CorrectionSurfaces(false);
|
||||
scaling_settings.StillsModulation(stills_modulation_flag);
|
||||
scaling_settings.StillsPartialityRefine(!simple_stills_flag);
|
||||
scaling_settings.ExpectedVarianceMerge(!no_expected_variance_merge);
|
||||
if (d_min_scale_merge)
|
||||
@@ -1446,7 +1410,6 @@ int main(int argc, char **argv) {
|
||||
// over-extrapolated under-captured fulls and, with the mosaicity fix, lifts rotation ISa/anomalous
|
||||
// substantially. Off for non-rot3d (no combine). An explicit --capture-uncertainty always wins.
|
||||
scaling_settings.CaptureUncertaintyCoeff(capture_uncertainty_arg.value_or(rotation_indexing ? 1.0 : 0.0));
|
||||
scaling_settings.PartialityUncertaintyCoeff(partiality_uncertainty_arg.value_or(0.0));
|
||||
scaling_settings.ForcedMosaicity(forced_mosaicity_arg);
|
||||
scaling_settings.MinCCForImage(min_image_cc / 100.0); // --min-image-cc is in percent; the setting is a fraction
|
||||
scaling_settings.OutlierRejectNsigma(
|
||||
|
||||
@@ -499,11 +499,6 @@ QWidget *JFJochViewerSettingsDock::BuildScalingSection() {
|
||||
corrections->setToolTip("Rotation only: fit a radiation-damage decay and a goniometer-frame absorption "
|
||||
"surface on the fulls. Cross-validated, so they no-op when their systematic is "
|
||||
"absent. On by default.");
|
||||
auto *modulation = new QCheckBox("Detector-plane modulation (stills)", this);
|
||||
modulation->setChecked(scaling_.GetStillsModulation());
|
||||
modulation->setToolTip("Stills: fit a detector-plane modulation (flat-field) surface over where each "
|
||||
"reflection lands, cross-validated so it no-ops when the systematic is absent. "
|
||||
"For rotation, modulation is part of \"Correction surfaces\" above.");
|
||||
auto *partRefine = new QCheckBox("Partiality post-refinement (stills)", this);
|
||||
partRefine->setChecked(scaling_.GetStillsPartialityRefine());
|
||||
partRefine->setToolTip("Stills: refine a per-crystal orientation tilt against the running merge and "
|
||||
@@ -515,23 +510,9 @@ QWidget *JFJochViewerSettingsDock::BuildScalingSection() {
|
||||
auto *highRes = new NumberLineEdit(0.3f, 5.0f, scaling_.GetHighResolutionLimit_A().value_or(2.0), 1, "Å", this);
|
||||
highRes->setEnabled(limitRes->isChecked());
|
||||
|
||||
// Stills partiality-uncertainty merge term: adds a systematic sigma ~c*(1-partiality)*<I> on partials,
|
||||
// so strong low-partiality partials are not over-trusted. Relevant once partials exist (partiality
|
||||
// post-refinement on, the default); the library auto-gates it to strong/medium data. 0 = off; ~2.5 rec.
|
||||
const double part_unc = scaling_.GetPartialityUncertaintyCoeff();
|
||||
auto *partUncertain = new QCheckBox("Partiality uncertainty", this);
|
||||
partUncertain->setChecked(part_unc > 0.0);
|
||||
partUncertain->setToolTip("Stills: add a systematic merge σ ~c·(1−partiality)·⟨I⟩ to partials so strong "
|
||||
"low-partiality partials are not over-trusted. Relevant with partiality "
|
||||
"post-refinement (the default); auto-gated to strong/medium data. ~2.5 "
|
||||
"recommended; unchecked = off.");
|
||||
auto *partUncertainCoeff = new NumberLineEdit(0.1f, 10.0f, part_unc > 0.0 ? part_unc : 2.5, 1, "", this);
|
||||
partUncertainCoeff->setEnabled(partUncertain->isChecked());
|
||||
|
||||
form->addRow("", friedel);
|
||||
form->addRow("", refineB);
|
||||
form->addRow("", corrections);
|
||||
form->addRow("", modulation);
|
||||
form->addRow("", partRefine);
|
||||
// Compact, and aligned with the checkboxes above: the limit checkbox + value sit together in the
|
||||
// field column (not as a row label, which would indent it differently).
|
||||
@@ -539,10 +520,6 @@ QWidget *JFJochViewerSettingsDock::BuildScalingSection() {
|
||||
resRow->addWidget(limitRes);
|
||||
resRow->addWidget(highRes, 1);
|
||||
form->addRow("", resRow);
|
||||
auto *puRow = new QHBoxLayout();
|
||||
puRow->addWidget(partUncertain);
|
||||
puRow->addWidget(partUncertainCoeff, 1);
|
||||
form->addRow("", puRow);
|
||||
section->setContentLayout(form);
|
||||
section->setExpanded(false); // folded on start (only geometry + unit cell start open)
|
||||
|
||||
@@ -550,24 +527,18 @@ QWidget *JFJochViewerSettingsDock::BuildScalingSection() {
|
||||
scaling_.MergeFriedel(friedel->isChecked());
|
||||
scaling_.RefineB(refineB->isChecked());
|
||||
scaling_.CorrectionSurfaces(corrections->isChecked());
|
||||
scaling_.StillsModulation(modulation->isChecked());
|
||||
scaling_.StillsPartialityRefine(partRefine->isChecked());
|
||||
scaling_.HighResolutionLimit_A(limitRes->isChecked()
|
||||
? std::optional<double>(highRes->value()) : std::nullopt);
|
||||
scaling_.PartialityUncertaintyCoeff(partUncertain->isChecked() ? partUncertainCoeff->value() : 0.0);
|
||||
emit scalingChanged(scaling_);
|
||||
};
|
||||
connect(friedel, &QCheckBox::toggled, this, [emitScaling] { emitScaling(); });
|
||||
connect(refineB, &QCheckBox::toggled, this, [emitScaling] { emitScaling(); });
|
||||
connect(corrections, &QCheckBox::toggled, this, [emitScaling] { emitScaling(); });
|
||||
connect(modulation, &QCheckBox::toggled, this, [emitScaling] { emitScaling(); });
|
||||
connect(partRefine, &QCheckBox::toggled, this, [emitScaling] { emitScaling(); });
|
||||
connect(limitRes, &QCheckBox::toggled, this, [emitScaling, highRes](bool on) {
|
||||
highRes->setEnabled(on); emitScaling(); });
|
||||
connect(highRes, &NumberLineEdit::newValue, this, [emitScaling] { emitScaling(); });
|
||||
connect(partUncertain, &QCheckBox::toggled, this, [emitScaling, partUncertainCoeff](bool on) {
|
||||
partUncertainCoeff->setEnabled(on); emitScaling(); });
|
||||
connect(partUncertainCoeff, &NumberLineEdit::newValue, this, [emitScaling] { emitScaling(); });
|
||||
return section;
|
||||
}
|
||||
|
||||
|
||||
Reference in New Issue
Block a user