Files
Jungfraujoch/common/Definitions.h
T
leonarski_fandClaude Opus 5 a6be35ccdb Azimuthal integration: optional sigma clipping of the reported profile
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>
2026-08-07 00:18:52 +02:00

73 lines
3.0 KiB
C++

// SPDX-FileCopyrightText: 2024 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
// SPDX-License-Identifier: GPL-3.0-only
#pragma once
#include "../fpga/pcie_driver/jfjoch_fpga.h"
#include <cstdint>
#include <cstddef>
#include <chrono>
#include <array>
#include <cmath>
constexpr float WVL_1A_IN_KEV = 12.39854f;
constexpr size_t CONVERTED_MODULE_LINES = 514;
constexpr size_t CONVERTED_MODULE_COLS = 1030;
constexpr size_t CONVERTED_MODULE_SIZE = CONVERTED_MODULE_LINES * CONVERTED_MODULE_COLS;
constexpr size_t JUNGFRAU_PACKET_SIZE_BYTES = 8192;
constexpr int MAX_IMAGE_NUMBER = 2*1024*1024;
constexpr std::chrono::nanoseconds MIN_COUNT_TIME = std::chrono::microseconds(3);
constexpr std::chrono::nanoseconds MIN_STORAGE_CELL_DELAY = std::chrono::nanoseconds(2100);
constexpr std::chrono::nanoseconds MIN_FRAME_TIME_JUNGFRAU_HALF_SPEED = std::chrono::microseconds(1000);
constexpr std::chrono::nanoseconds MIN_FRAME_TIME_JUNGFRAU_FULL_SPEED = std::chrono::microseconds(470);
constexpr std::chrono::nanoseconds MIN_FRAME_TIME_EIGER = std::chrono::microseconds(250);
constexpr std::chrono::nanoseconds MAX_COUNT_TIME_JUNGFRAU = std::chrono::microseconds(2000);
constexpr std::chrono::nanoseconds FRAME_TIME_PEDE_G1G2 = std::chrono::microseconds(10*1000);
constexpr std::chrono::nanoseconds PSI_JUNGFRAU_READOUT_TIME = std::chrono::microseconds(20);
constexpr std::chrono::nanoseconds PSI_EIGER_READOUT_TIME = std::chrono::microseconds(20);
constexpr std::chrono::nanoseconds DARK_MASK_FRAME_TIME = std::chrono::milliseconds(10);
constexpr float MIN_ENERGY_KEV = 0.001;
constexpr float MAX_ENERGY_KEV = 500.0;
constexpr float DEFAULT_G0_FACTOR = 41.0f;
constexpr float DEFAULT_G1_FACTOR = -1.439f;
constexpr float DEFAULT_G2_FACTOR = -0.1145f;
constexpr float DEFAULT_HG0_FACTOR = 100.0f;
constexpr int MAX_SPOT_COUNT = 64 * 1024;
constexpr uint32_t MASK_PEDESTAL_G0_RMS_LIMIT = (1U<<4);
constexpr size_t PEDESTAL_MIN_IMAGE_COUNT = 128;
constexpr uint16_t PEDESTAL_WRONG = (UINT16_MAX);
constexpr size_t PEDESTAL_G0_WRONG_GAIN_ALLOWED_COUNT = 2;
constexpr size_t MESSAGE_SIZE_FOR_START_END = (256*1024*1024); // pessimistic highest value
constexpr float LAB6_CELL_A = 4.156468f;
// Ice ring resolution taken from:
// Moreau, Atakisi, Thorne, Acta Cryst D77, 2021, 540,554
// https://journals.iucr.org/d/issues/2021/04/00/tz5104/index.html
constexpr std::array<float, 11> ICE_RING_RES_A = {3.895, 3.661, 3.438, 2.667, 2.249, 2.068, 1.947, 1.916, 1.882, 1.719, 1.522};
// True when resolution d (Angstrom) sits within half_width of a hexagonal-ice powder ring, in the
// q = 2*pi/d units the spot-finder uses (ice_ring_width_Q_recipA). Used to drop ice-contaminated
// reflections from scaling/merging when ice-ring handling is enabled.
inline bool IsOnIceRing(float d_A, float half_width_q_recipA) {
if (!(d_A > 0.0f))
return false;
constexpr float two_pi = 6.283185307f;
const float q = two_pi / d_A;
for (const float ice_d : ICE_RING_RES_A)
if (std::fabs(q - two_pi / ice_d) < half_width_q_recipA)
return true;
return false;
}