ice score: measure against a smooth whole-profile background, not shoulders
The score's baseline was two adjacent shoulder bins with a bin-overlap bug - the ring's edge bins were counted in both the ring and the shoulder, since GetMeanValueOfBins is inclusive. At the typical (coarse) azint binning (dq ~ 0.05 in q, wider than the 0.03 ring half-width) a shoulder is only ~1 bin, so the ratio was noisy and poorly separated. Replace it with the ring intensity over a smooth whole-profile background: a running median of the non-ice bins, interpolated under each ring. Clean crystals now sit at ~1.0 and ice separates far more cleanly on /data/rotation_test: cytC 1.06->1.03, lysoC 1.23->2.77, EP_cs_01-17 1.67->4.51 (max 11.4). A z-score / abnormality probability was tried but is uninformative here - with many photons any real ice ring is highly significant, so the useful discriminator is the ice magnitude (this ratio), noted in CPU_DATA_ANALYSIS.md. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
This commit is contained in:
Claude Fable 5
@@ -125,7 +125,7 @@ public:
|
||||
bool spotsIceIsSet() const;
|
||||
void unsetSpots_ice();
|
||||
/// <summary>
|
||||
/// Strongest hexagonal-ice ring band/shoulder intensity ratio (1 = no ice)
|
||||
/// Strongest hexagonal-ice ring intensity over the smooth radial background (1 = no ice)
|
||||
/// </summary>
|
||||
float getIceRingScore() const;
|
||||
void setIceRingScore(float const value);
|
||||
|
||||
@@ -1582,7 +1582,7 @@ components:
|
||||
ice_ring_score:
|
||||
type: number
|
||||
format: float
|
||||
description: Strongest hexagonal-ice ring band/shoulder intensity ratio (1 = no ice)
|
||||
description: Strongest hexagonal-ice ring intensity over the smooth radial background (1 = no ice)
|
||||
index:
|
||||
type: integer
|
||||
format: int64
|
||||
|
||||
@@ -176,23 +176,74 @@ float AzimuthalIntegrationProfile::GetBkgEstimate(const AzimuthalIntegrationSett
|
||||
|
||||
float AzimuthalIntegrationProfile::GetIceRingScore(const AzimuthalIntegrationSettings &settings,
|
||||
float half_width_q) const {
|
||||
// For each hexagonal-ice powder ring, the mean profile intensity in the ring band (+/- half_width in
|
||||
// q = 2*pi/d) over a baseline from the two shoulders just outside it. Report the strongest ring's
|
||||
// ratio (1 = no ice); rings whose band+shoulders fall off the measured q-range are skipped.
|
||||
// Strongest hexagonal-ice ring's intensity relative to the background *under* it (1 = no ice). The
|
||||
// background is a smooth whole-profile estimate: a running median of the NON-ice bins, interpolated to
|
||||
// each ring position - not a couple of adjacent shoulder bins (the azint binning is coarser than the
|
||||
// ring width, so a local shoulder is only ~1 bin and a narrow ratio is noisy and can double-count the
|
||||
// ring's own edge). Clean profiles then sit at ~1 at every ring; ice makes the ring bin stand out.
|
||||
constexpr float two_pi = 6.283185307f;
|
||||
const std::vector<float> prof = GetResult1D();
|
||||
const int nq = static_cast<int>(prof.size());
|
||||
const float low_q = settings.GetLowQ_recipA();
|
||||
const float high_q = settings.GetHighQ_recipA();
|
||||
const float dq = settings.GetQSpacing_recipA();
|
||||
if (nq < 12 || !(dq > 0.0f))
|
||||
return 1.0f;
|
||||
|
||||
auto q_of = [&](int i) { return low_q + (static_cast<float>(i) + 0.5f) * dq; };
|
||||
auto on_ice = [&](float q) {
|
||||
for (const float d : ICE_RING_RES_A)
|
||||
if (std::fabs(q - two_pi / d) < 1.5f * half_width_q)
|
||||
return true;
|
||||
return false;
|
||||
};
|
||||
|
||||
// Non-ice, finite, positive bins (ascending q) carry the background.
|
||||
std::vector<int> base;
|
||||
for (int i = 0; i < nq; ++i)
|
||||
if (std::isfinite(prof[i]) && prof[i] > 0.0f && !on_ice(q_of(i)))
|
||||
base.push_back(i);
|
||||
if (base.size() < 8)
|
||||
return 1.0f;
|
||||
|
||||
// Running median over the base bins => a smooth background robust to the ice peaks.
|
||||
constexpr int K = 4;
|
||||
std::vector<float> base_bg(base.size());
|
||||
std::vector<float> window;
|
||||
for (int j = 0; j < static_cast<int>(base.size()); ++j) {
|
||||
const int lo = std::max(0, j - K);
|
||||
const int hi = std::min(static_cast<int>(base.size()), j + K + 1);
|
||||
window.clear();
|
||||
for (int m = lo; m < hi; ++m)
|
||||
window.push_back(prof[base[m]]);
|
||||
std::sort(window.begin(), window.end());
|
||||
base_bg[j] = window[window.size() / 2];
|
||||
}
|
||||
|
||||
float score = 1.0f;
|
||||
for (const float ice_d : ICE_RING_RES_A) {
|
||||
const float q = two_pi / ice_d;
|
||||
if (q - 2 * half_width_q < low_q || q + 2 * half_width_q > high_q)
|
||||
for (const float d : ICE_RING_RES_A) {
|
||||
const float qr = two_pi / d;
|
||||
const int b = static_cast<int>(std::lround((qr - low_q) / dq - 0.5f));
|
||||
if (b < 0 || b >= nq || !std::isfinite(prof[b]) || prof[b] <= 0.0f)
|
||||
continue;
|
||||
const float ring = GetMeanValueOfBins(settings.QToBin(q - half_width_q), settings.QToBin(q + half_width_q));
|
||||
const float lo = GetMeanValueOfBins(settings.QToBin(q - 2 * half_width_q), settings.QToBin(q - half_width_q));
|
||||
const float hi = GetMeanValueOfBins(settings.QToBin(q + half_width_q), settings.QToBin(q + 2 * half_width_q));
|
||||
const float baseline = 0.5f * (lo + hi);
|
||||
if (std::isfinite(baseline) && baseline > 0.0f && std::isfinite(ring))
|
||||
score = std::max(score, ring / baseline);
|
||||
// Linear-interpolate the smooth background to the ring position.
|
||||
float bg;
|
||||
if (qr <= q_of(base.front()))
|
||||
bg = base_bg.front();
|
||||
else if (qr >= q_of(base.back()))
|
||||
bg = base_bg.back();
|
||||
else {
|
||||
bg = NAN;
|
||||
for (int j = 0; j + 1 < static_cast<int>(base.size()); ++j) {
|
||||
const float qa = q_of(base[j]), qb = q_of(base[j + 1]);
|
||||
if (qa <= qr && qr <= qb) {
|
||||
const float t = (qb > qa) ? (qr - qa) / (qb - qa) : 0.0f;
|
||||
bg = base_bg[j] + t * (base_bg[j + 1] - base_bg[j]);
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
if (std::isfinite(bg) && bg > 0.0f)
|
||||
score = std::max(score, prof[b] / bg);
|
||||
}
|
||||
return score;
|
||||
}
|
||||
|
||||
@@ -44,8 +44,8 @@ public:
|
||||
|
||||
float GetMeanValueOfBins(uint16_t min_bin, uint16_t max_bin) const;
|
||||
float GetBkgEstimate(const AzimuthalIntegrationSettings& settings) const;
|
||||
// Single per-image ice indicator: the strongest hexagonal-ice ring's band/shoulder intensity ratio
|
||||
// (1 = no ice, >1 = ice above the local background). Max over the rings in range; 1 if none.
|
||||
// Single per-image ice indicator: the strongest hexagonal-ice ring's intensity relative to the smooth
|
||||
// radial background interpolated under it (1 = no ice, >1 = ice). See the .cpp for the background fit.
|
||||
float GetIceRingScore(const AzimuthalIntegrationSettings& settings, float half_width_q) const;
|
||||
MultiLinePlot GetPlot(bool force_1d = false, PlotAzintUnit plot_unit = PlotAzintUnit::Q_recipA) const;
|
||||
AzimuthalIntegrationProfile& operator+=(const AzimuthalIntegrationProfile& profile); // Not thread safe
|
||||
|
||||
@@ -119,7 +119,7 @@ struct DataMessage {
|
||||
std::vector<uint64_t> az_int_profile_count;
|
||||
|
||||
std::optional<float> bkg_estimate;
|
||||
std::optional<float> ice_ring_score; // strongest ice-ring band/shoulder ratio from the azint profile (1 = none)
|
||||
std::optional<float> ice_ring_score; // strongest ice ring over the smooth radial background (1 = none)
|
||||
|
||||
std::optional<bool> indexing_result;
|
||||
std::optional<CrystalLattice> indexing_lattice;
|
||||
|
||||
+2
-2
@@ -215,7 +215,7 @@ See [DECTRIS documentation](https://github.com/dectris/documentation/tree/main/s
|
||||
| packets_expected | uint64 | Number of packets expected per image (in units of 2 kB) | | |
|
||||
| packets_received | uint64 | Number of packets received per image (in units of 2 kB) | | |
|
||||
| bkg_estimate | float | Mean value for pixels in resolution range from 3.0 to 5.0 A \[photons\] | | |
|
||||
| ice_ring_score | float | Strongest hexagonal-ice ring band/shoulder intensity ratio from the azint profile (1 = no ice) | | |
|
||||
| ice_ring_score | float | Strongest hexagonal-ice ring intensity over the smooth radial background (1 = no ice) | | |
|
||||
| beam_corr_x | float | Beam center correction X applied during processing \[pixel\] | | X |
|
||||
| beam_corr_y | float | Beam center correction Y applied during processing \[pixel\] | | X |
|
||||
| image_scale_factor | float | Scaling result: Image scale factor (g) | | X |
|
||||
@@ -306,7 +306,7 @@ See [DECTRIS documentation](https://github.com/dectris/documentation/tree/main/s
|
||||
| spot_count_indexed | Array(int32) | Per-image number of spots fitting indexing solution | |
|
||||
| image_indexed | Array(uint8) | Per-image indexing result; 0 = not indexed, nonzero = indexed | |
|
||||
| v_bkg_estimate | Array(float) | Per-image background estimate | |
|
||||
| ice_ring_score | Array(float) | Per-image strongest ice-ring band/shoulder intensity ratio (1 = no ice) | |
|
||||
| ice_ring_score | Array(float) | Per-image strongest ice-ring intensity over the smooth radial background (1 = no ice) | |
|
||||
| profile_radius | Array(float) | Per-image profile radius \[Angstrom^-1\] | |
|
||||
| mosaicity | Array(float) | Per-image mosaicity \[degree\] | |
|
||||
| bFactor | Array(float) | Per-image estimated B-factor \[Angstrom^2\] | |
|
||||
|
||||
@@ -180,7 +180,7 @@ Special cases:
|
||||
|
||||
Spot finding can be restricted to a resolution range $[d_\mathrm{high}, d_\mathrm{low}]$ by masking pixels outside the range. Optionally, pixels in identified ice-ring regions can be tagged so that subsequent indexing/refinement may include or exclude them (see §4 and §6).
|
||||
|
||||
A single per-image **ice-ring score** is derived from the azimuthally-integrated radial profile: for each hexagonal-ice powder ring (positions $d$ from Moreau *et al.*, Acta Cryst D77, 2021), the mean profile intensity in the ring band ($\pm$ ice-ring half-width in $q$) is divided by a baseline interpolated from the two shoulders just outside it, and the strongest ring's ratio is reported (1 = no ice, $>1$ = ice above background). It is stored per image (`ice_ring_score`, HDF5 `/entry/MX/iceRingScore`) as a monitoring quantity. Note this is distinct from the merge-time ice masking, which is data-driven from the per-ring merged CC1/2 rather than this background ratio.
|
||||
A single per-image **ice-ring score** is derived from the azimuthally-integrated radial profile: for each hexagonal-ice powder ring (positions $d$ from Moreau *et al.*, Acta Cryst D77, 2021), the profile intensity at the ring is divided by a smooth background estimated from the *whole* profile — a running median of the non-ice bins, interpolated under each ring — and the strongest ring's ratio is reported (1 = no ice, $>1$ = ice above background). A whole-profile background is used rather than a couple of adjacent shoulder bins because the radial binning is typically coarser than the ring width (a local shoulder is only ~1 bin), so a narrow ratio is noisy and can double-count the ring's own edge. (A significance/z-score was considered but is uninformative here: with many photons any real ice ring is highly significant, so the discriminating quantity is the ice *magnitude*, i.e. this ratio.) It is stored per image (`ice_ring_score`, HDF5 `/entry/MX/iceRingScore`) as a monitoring quantity, distinct from the merge-time ice masking, which is data-driven from the per-ring merged CC1/2.
|
||||
|
||||
A further optional safeguard removes isolated high-resolution “spur” spots by detecting large gaps in $1/d$ (or $q$) space and discarding spots beyond the gap. This is intended for macromolecular diffraction where edge-of-detector backgrounds can be extremely low.
|
||||
|
||||
|
||||
+1
-1
@@ -251,7 +251,7 @@ In legacy/VDS mode these live in the data files and are linked/virtual-stacked i
|
||||
| `resolutionEstimate` | Å | diffraction resolution estimate |
|
||||
| `integratedReflections` | | number of integrated reflections |
|
||||
| `bkgEstimate` | photons | mean background in the 3–5 Å resolution band |
|
||||
| `iceRingScore` | ratio | strongest hexagonal-ice ring band/shoulder intensity ratio (1 = no ice) |
|
||||
| `iceRingScore` | ratio | strongest hexagonal-ice ring intensity over the smooth radial background (1 = no ice) |
|
||||
| `beam_corr_x`, `beam_corr_y` | pixel | beam-center correction applied during processing |
|
||||
| `imageScaleFactor` | | on-the-fly per-image scale factor *g* |
|
||||
| `imageScaleCC` | | on-the-fly scaling correlation coefficient |
|
||||
|
||||
@@ -1011,7 +1011,7 @@ export type scan_result = {
|
||||
*/
|
||||
spots_ice?: number;
|
||||
/**
|
||||
* Strongest hexagonal-ice ring band/shoulder intensity ratio (1 = no ice)
|
||||
* Strongest hexagonal-ice ring intensity over the smooth radial background (1 = no ice)
|
||||
*/
|
||||
ice_ring_score?: number;
|
||||
/**
|
||||
|
||||
Reference in New Issue
Block a user