Files
Jungfraujoch/common/Definitions.h
T
leonarski_fandClaude Opus 5 6516bc96af Ice rings: carry the list past 1.5 A, where ice does not stop
The eleven measured bands end at 1.522 A because their source says so in its
own words - "pure hexagonal ice has 11 diffraction rings between 4 and 1.5 A
resolution" - and its subject was detecting ice in deposited data, not masking
it. On a detector that reaches further, the rings it does not list are the ones
left in the data: on the strong rotation set just added to the battery, 44% of
every image's spots sit in ice bands, and beyond 1.5 A the spot list is ice and
nothing else, which is why the resolution estimate read the ice rather than the
crystal.

There is nothing measured to copy below 1.522 A, so the eight added bands are
calculated. Enumerating hkl is not enough and the code already said so: ice Ih
is P6_3/mmc with O on 4f, and most of what enumeration emits is extinguished by
the OXYGEN SUBLATTICE rather than by the space group - which is why (004) at
1.830 A and (104) at 1.657 A are missing from the measured list although they
sit inside its range and its reflection conditions allow them (for (00l) the
structure factor goes as cos(2*pi*l*z), and z ~ 1/16 kills l = 4). So compute
structure factors - oxygen only, the hydrogens being half-occupancy disordered
and weak to X-rays - and keep the lines reaching 3% of the strongest. That rule
REPRODUCES THE MEASURED ELEVEN EXACTLY and every line it drops inside their
range computes to zero, which is what makes it trustworthy below 1.522 A. It
stops at 1.170 A: below that the real lines fall to 2-3% while the extinct ones
rise to about 1%, and an oxygen-only calculation cannot separate them honestly.

Every added band was independently confirmed in the data - the spot-count
histogram of the strong set peaks at each of them and is empty between - and
every line the rule calls extinct is absent there too.

Costs, measured. The bands are inert above 1.6 A: on 38 of 39 battery sets the
profile ice score does not move at all, and the one that appeared to (a jet set,
1.25 -> 2.60) does not on the peak-excluded profile the score actually uses -
that was Bragg peaks in the plain profile, which is what the peak exclusion is
for. Where a detector does reach past 1.5 A the bands cover more of reciprocal
space: unchanged at 1.6 A, +7.4 points at 1.4 A, +16.5 at 1.18 A. On the strong
set that is 17% -> 27% of reflections held out of the scale fit, and it shows:
the spot resolution estimate improves from 1.33 to 1.46 A against a truth near
1.42, while CC1/2 falls 98.5 -> 97.6% and ISa 3.58 -> 3.37. Ice handling only
runs at all on a run that trips the ice gate, so a clean crystal pays nothing.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01FBumeJVx4oeXxiBRpkrE5H
2026-08-28 11:27:28 +02:00

99 lines
5.1 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;
// Defaults for images_per_file when the acquisition does not ask for a particular number; see
// DiffractionExperiment::GetImagesPerFile.
constexpr int64_t DEFAULT_IMAGES_PER_FILE = 1000;
// Above this, a rotation sweep is split rather than written to one file.
constexpr int64_t ROTATION_SINGLE_FILE_IMAGE_LIMIT = 20000;
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;
// Hexagonal-ice ring positions. The first eleven are the measured ones from:
// Moreau, Atakisi, Thorne, Acta Cryst D77, 2021, 540,554
// https://journals.iucr.org/d/issues/2021/04/00/tz5104/index.html
// That paper's Table 1 stops at 1.522 A because it says so in its own words - "pure hexagonal ice has
// 11 diffraction rings between 4 and 1.5 A resolution" - and its subject was detecting ice in the PDB,
// not masking it. Ice does not stop there, and on a detector that reaches past 1.5 A the rings it does
// not list are the ones left in the data: one strongly diffracting set had 44% of every frame's spots
// in ice bands, and beyond 1.5 A its spots were ice and nothing else.
//
// The rest are calculated, because past that paper there is nothing measured to copy. Ice Ih is
// P6_3/mmc with O on 4f, so enumerating hkl is not enough - most of what it emits is extinguished by
// the oxygen sublattice rather than by the space group, which is why (004) at 1.830 A and (104) at
// 1.657 A are absent from the measured list although they lie inside its range and its reflection
// conditions allow them. Structure factors were computed instead (oxygen only - the hydrogens are
// half-occupancy disordered and scatter X-rays weakly) and the lines kept are those reaching 3% of the
// strongest. That reproduces the measured eleven exactly, and every line it drops in their range
// computes to zero, which is what makes the same rule trustworthy below 1.522 A.
// Following Roettger, Endriss, Ihringer, Doyle & Kuhs (1994) Acta Cryst. B50, 644-648 for the cell.
//
// It stops at 1.170 A: below that the calculated real lines fall to 2-3% while the extinct ones rise
// to about 1%, and an oxygen-only calculation cannot separate them honestly any further.
constexpr std::array<float, 19> 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,
1.472, 1.443, 1.371, 1.366, 1.298, 1.261, 1.224,
1.170};
// 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;
}