The profile is the MEAN of each bin, so a few strong reflections landing in a bin lift it exactly as a smooth powder ring does. That is the wrong quantity whenever the profile is wanted as a background rather than as a measurement of what is in the bin - the ice score being the case in point, where reading a plain profile INVERTED the metric: over 37 rotation crystals the two highest-scoring crystals had no ice at all. The adaptive spot finder already computes the right thing, a sigma-clipped per-resolution-ring background, as a byproduct of its own threshold. Where it runs, the ice score uses that. Where it does not - --no-adaptive-spots, --azint-only, and anything reading the profile the broker wrote - there was no way to get it. This adds one: azim_int_settings.sigma_clip (rugnux --azim-sigma-clip), 0 = off, minimum 2 because a tighter clip rejects a large part of a clean Gaussian bin and biases the estimate low rather than removing outliers. Two clip passes follow the plain one, matching the finder's recipe - the first pass's standard deviation is itself inflated by the peaks being removed, so one pass leaves a threshold that is still too generous. A bin with fewer than eight pixels is left alone: at the detector edge and behind the beam stop there is no spread to clip on. Both engines do it. On the GPU the accept range is computed by a small kernel and stays resident, so a clip pass is one more read of the same pixels and no round trip; the two accumulation kernels take the range as a pointer that is null on the plain pass. Measured on a JUNGFRAU rotation dataset, non-adaptive path: azimuthal integration 0.02 -> 0.06 ms per image, exactly the 3x the extra passes predict, against a 0.34 ms per-image total. Note what the result IS: the smooth background under the peaks, not the bin mean. It should not be switched on where a ring's integrated intensity is wanted - the powder-ring geometry fit reads ring peaks, and those are what a clip is designed to remove. Off by default, so nothing changes unless it is asked for. Not exposed over the REST API - that needs the generated model regenerated, which is a separate step. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
295 lines
8.7 KiB
C++
295 lines
8.7 KiB
C++
// SPDX-FileCopyrightText: 2026 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
|
|
// SPDX-License-Identifier: GPL-3.0-only
|
|
|
|
|
|
#include "ScalingSettings.h"
|
|
|
|
|
|
ScalingSettings& ScalingSettings::MergeFriedel(bool input) {
|
|
merge_friedel = input;
|
|
return *this;
|
|
}
|
|
|
|
ScalingSettings& ScalingSettings::HighResolutionLimit_A(double limit) {
|
|
if (limit <= 0.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterBelowMin, "High resolution limit must be positive");
|
|
high_resolution_limit_A = limit;
|
|
return *this;
|
|
}
|
|
|
|
ScalingSettings& ScalingSettings::HighResolutionLimit_A(std::optional<double> limit) {
|
|
if (limit.has_value() && limit.value() <= 0.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterBelowMin, "High resolution limit must be positive");
|
|
high_resolution_limit_A = limit;
|
|
return *this;
|
|
}
|
|
|
|
bool ScalingSettings::GetMergeFriedel() const {
|
|
return merge_friedel;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::RefineRotationWedge(bool input) {
|
|
refine_wedge = input;
|
|
return *this;
|
|
}
|
|
|
|
bool ScalingSettings::GetRefineWedge() const {
|
|
return refine_wedge;
|
|
}
|
|
|
|
std::optional<double> ScalingSettings::GetHighResolutionLimit_A() const {
|
|
return high_resolution_limit_A;
|
|
}
|
|
|
|
double ScalingSettings::GetMinMosaicity() const {
|
|
return 0.001;
|
|
}
|
|
|
|
double ScalingSettings::GetMaxMosaicity() const {
|
|
return 1.0;
|
|
}
|
|
|
|
double ScalingSettings::GetMinWedge() const {
|
|
return 0.001;
|
|
|
|
}
|
|
|
|
double ScalingSettings::GetMaxWedge() const {
|
|
return 10.0;
|
|
}
|
|
|
|
double ScalingSettings::GetDefaultMosaicity() const {
|
|
return 0.1;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::RotationWedgeForScaling(std::optional<double> input) {
|
|
if (input) {
|
|
// TODO: Use fmt
|
|
if (input.value() < GetMinWedge() || input.value() > GetMaxWedge())
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
|
"Wedge for scaling must be between " + std::to_string(GetMinWedge()) +
|
|
" and " + std::to_string(GetMaxWedge()));
|
|
}
|
|
wedge_for_scaling = input;
|
|
return *this;
|
|
}
|
|
|
|
std::optional<double> ScalingSettings::GetRotationWedgeForScaling() const {
|
|
return wedge_for_scaling;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::MinPartiality(double input) {
|
|
if (min_partiality < 0.0 || min_partiality > 1.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Min partiality must be between 0 and 1");
|
|
min_partiality = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetMinCCForImage() const {
|
|
return min_cc_for_image;
|
|
}
|
|
|
|
double ScalingSettings::GetSearchMinZeta() const {
|
|
return search_min_zeta;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::SearchMinZeta(double input) {
|
|
if (input < 0.0 || input >= 1.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
|
"Search zeta limit must be in [0,1)");
|
|
search_min_zeta = input;
|
|
return *this;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::MinCCForImage(double input) {
|
|
if (input < 0.0 || input > 1.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Min CC for image must be between 0 and 1");
|
|
min_cc_for_image = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetOutlierRejectNsigma() const {
|
|
return outlier_reject_nsigma;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::OutlierRejectNsigma(double input) {
|
|
outlier_reject_nsigma = input; // <= 0 disables; no upper bound (large = effectively off)
|
|
return *this;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::ScaleFulls(bool input) {
|
|
scale_fulls = input;
|
|
return *this;
|
|
}
|
|
|
|
bool ScalingSettings::GetScaleFulls() const {
|
|
return scale_fulls;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::AbsorptionIter(int input) {
|
|
if (input < 0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Absorption iterations must be non-negative");
|
|
absorption_iter = input;
|
|
return *this;
|
|
}
|
|
|
|
int ScalingSettings::GetAbsorptionIter() const {
|
|
return absorption_iter;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::CorrectionSurfaces(bool input) {
|
|
correction_surfaces = input;
|
|
return *this;
|
|
}
|
|
|
|
bool ScalingSettings::GetCorrectionSurfaces() const {
|
|
return correction_surfaces;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::StillsPartialityRefine(bool input) {
|
|
stills_partiality_refine = input;
|
|
return *this;
|
|
}
|
|
|
|
bool ScalingSettings::GetStillsPartialityRefine() const {
|
|
return stills_partiality_refine;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::ExpectedVarianceMerge(bool input) {
|
|
expected_variance_merge = input;
|
|
return *this;
|
|
}
|
|
|
|
bool ScalingSettings::GetExpectedVarianceMerge() const {
|
|
return expected_variance_merge;
|
|
}
|
|
|
|
float ScalingSettings::GetIceMinScore() const {
|
|
return ice_min_score;
|
|
}
|
|
|
|
float ScalingSettings::GetIceMinSpotRatio() const {
|
|
return ice_min_spot_ratio;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::IceMinSpotRatio(float input) {
|
|
if (input < 0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Ice spot-ratio gate must be non-negative");
|
|
ice_min_spot_ratio = input;
|
|
return *this;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::IceMinScore(float input) {
|
|
if (input < 0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Ice score gate must be non-negative");
|
|
ice_min_score = input;
|
|
return *this;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::SmoothGDegrees(double input) {
|
|
if (input < 0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Smooth-G range must be non-negative");
|
|
smooth_g_deg = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetSmoothGDegrees() const {
|
|
return smooth_g_deg;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::RelativeBDegrees(double input) {
|
|
if (input < 0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "Relative-B batch width must be non-negative");
|
|
relative_b_deg = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetRelativeBDegrees() const {
|
|
return relative_b_deg;
|
|
}
|
|
|
|
double ScalingSettings::GetMinPartiality() const {
|
|
return min_partiality;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::ForcedMosaicity(std::optional<double> input) {
|
|
if (input.has_value() && (input.value() < GetMinMosaicity() || input.value() > GetMaxMosaicity()))
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
|
"Forced mosaicity must be between " + std::to_string(GetMinMosaicity()) +
|
|
" and " + std::to_string(GetMaxMosaicity()));
|
|
forced_mosaicity = input;
|
|
return *this;
|
|
}
|
|
|
|
std::optional<double> ScalingSettings::GetForcedMosaicity() const {
|
|
return forced_mosaicity;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::CaptureUncertaintyCoeff(double input) {
|
|
if (input < 0.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
|
"Capture uncertainty coefficient must be non-negative");
|
|
capture_uncertainty_coeff = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetCaptureUncertaintyCoeff() const {
|
|
return capture_uncertainty_coeff;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::MinCapturedFraction(double input) {
|
|
if (input < 0.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
|
"Minimum captured fraction must be non-negative");
|
|
min_captured_fraction = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetMinCapturedFraction() const {
|
|
return min_captured_fraction;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::RfreeFraction(double input) {
|
|
if (input < 0.0 || input > 1.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, "R-free fraction must be between 0 and 1");
|
|
rfree_fraction = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetRfreeFraction() const {
|
|
return rfree_fraction;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::ResolutionCutoff(ResolutionCutoffMethod input) {
|
|
resolution_cutoff = input;
|
|
return *this;
|
|
}
|
|
|
|
ResolutionCutoffMethod ScalingSettings::GetResolutionCutoff() const {
|
|
return resolution_cutoff;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::ResolutionCCTarget(double input) {
|
|
if (input <= 0.0 || input >= 1.0)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
|
"Resolution CC target must be between 0 and 1");
|
|
resolution_cc_target = input;
|
|
return *this;
|
|
}
|
|
|
|
double ScalingSettings::GetResolutionCCTarget() const {
|
|
return resolution_cc_target;
|
|
}
|
|
|
|
ScalingSettings &ScalingSettings::ReportShellCount(int input) {
|
|
if (input < 1)
|
|
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
|
|
"Number of report shells must be at least 1");
|
|
report_shell_count = input;
|
|
return *this;
|
|
}
|
|
|
|
int ScalingSettings::GetReportShellCount() const {
|
|
return report_shell_count;
|
|
}
|