Files
Jungfraujoch/common/DatasetSettings.cpp
leonarski_f a39fd29f77
Build Packages / XDS test (JFJoch plugin) (push) Successful in 11m4s
Build Packages / Unit tests (push) Skipped
Build Packages / build:windows:nocuda (push) Successful in 17m46s
Build Packages / build:windows:cuda (push) Successful in 20m20s
Build Packages / build:viewer-tgz:cpu (push) Successful in 15m56s
Build Packages / build:viewer-tgz:cuda (push) Successful in 17m57s
Build Packages / build:rugnux-tgz (x86_64) (push) Successful in 14m10s
Build Packages / build:rugnux:windows (push) Successful in 11m12s
Build Packages / build:rugnux:aarch64 (cross) (push) Successful in 7m14s
Build Packages / build:rpm (rocky8_nocuda) (push) Successful in 22m13s
Build Packages / build:rpm (rocky9_nocuda) (push) Successful in 19m17s
Build Packages / build:rpm (ubuntu2204_nocuda) (push) Successful in 21m21s
Build Packages / build:rpm (ubuntu2404_nocuda) (push) Successful in 17m26s
Build Packages / build:rpm (rocky8_sls9) (push) Successful in 23m56s
Build Packages / build:rpm (rocky9_sls9) (push) Successful in 20m48s
Build Packages / build:rpm (rocky8) (push) Successful in 23m43s
Build Packages / build:rpm (rocky9) (push) Successful in 20m38s
Build Packages / build:rpm (ubuntu2204) (push) Successful in 24m57s
Build Packages / build:rpm (ubuntu2404) (push) Successful in 20m58s
Build Packages / XDS test (durin plugin) (push) Successful in 10m43s
Build Packages / Generate python client (push) Successful in 47s
Build Packages / Build documentation (push) Successful in 1m5s
Build Packages / Create release (push) Skipped
Build Packages / XDS test (neggia plugin) (push) Successful in 8m57s
Build Packages / DIALS test (push) Successful in 18m40s
v1.0.0-rc.167 (#77)
* `rugnux --model` reports CC(model, data) - the correlation of the merged intensities with the placed, scaled model - by resolution shell, on the same shells as CC1/2, with the reflection count and a significance for each.
* `rugnux --model` fits the model's scale, anisotropic B and bulk-solvent parameters on the working reflections only, so the R-free it reports is measured against a model no free reflection helped scale.
* The bulk-solvent parameters of `rugnux --model` are searched over their physically meaningful range instead of being fitted without bounds, so a model is never scaled with a solvent term that has silently switched itself off.
* The rigid-body placement of `rugnux --model` uses the same bounded bulk solvent as the reported fit, so a model is no longer placed against a target carrying a solvent term with no physical meaning.
* `rugnux --model` puts the model into the data's own description of the lattice before placing it, so a model whose cell is written on other axes - I-centred where the run indexed C-centred, a different unique axis, a permuted orthorhombic cell - is placed rather than scored where it was read; `MODEL_CHANGE_OF_BASIS=` and `MODEL_SETTING_AS_READ=` report it when it happens.
* The rugnux results report opens with a summary - `VERDICT=` (`OK`, `WARNINGS`, `UNUSABLE`, `FAILED`), `VERDICT_TEXT=`, `PATHOLOGY_FLAGS=` with one closed-vocabulary code per condition that warned, and the `WARNING:` lines, which used to close the file - and the sections after it are renumbered 1-5 with no gaps.
* `rugnux --developer` writes the full results report - the pipeline-internal keys and the long explanations the default report now leaves out - and `--finalist-ledger` adds the evidence for every space group the search considered, not only the one it adopted.
* The results report warns when the merged data carry no usable signal and when too little of reciprocal space was measured inside the fitted resolution, and omits `FITTED_RESOLUTION` where the CC1/2 curve it is fitted on never falls off.
* rugnux detects translational pseudo-symmetry and reports it under the `PSEUDO_TRANSLATION` flag as `TNCS_DETECTED=` and the `TNCS_*` keys - a translation the merged data are exactly invariant under is reported as `UNDECLARED_LATTICE_TRANSLATION=` under `LATTICE_TRANSLATION` instead - and a detected pseudo-translation can no longer buy a false screw axis in the space-group search or hide a twin from the L-test (`L_TEST_VS_TNCS=`).
* The space-group search determines glide planes from zonal systematic absences, so a non-Sohncke space group such as P 2_1/c or Pbca is named where the run previously stopped at its Sohncke subgroup; `SOHNCKE_SPACE_GROUP=` carries the best Sohncke group beside it on every run that searched, and a centre of symmetry is never claimed.
* Where the cell metric carries more rotational symmetry than the Bravais class the indexer named, the extra rotations are put to the intensities and the space-group search is asked again on the metric's own cell - adopted only where the intensities confirm the higher symmetry - so a lattice that is nearly but not exactly hexagonal, or whose reduction landed in a sub-cell, still reaches its true point group.
* Systematic-absence calls rest on the evidence rather than on counts: a screw axis whose absent class the data show extinct is no longer refused because a handful of reflections in it read as present, and `SPACE_GROUP_ALTERNATIVES=` no longer drops a candidate that differs only on a zone the sweep never measured.
* A reference correlation measured on too few reflections is refused instead of scored zero, so a run given a reference MTZ is no longer reindexed on an operator that mapped almost everything outside the reference's coverage.
* A frame counts as indexed from 6 spots on its lattice rather than 9, so a weakly diffracting crystal whose frames cannot carry 9 is no longer refused the lattice it fits; `--min-indexed-spots` overrides it.
* `-C` accepts a known cell in any equivalent description - conventional or primitive, centred or not - instead of only the reduced primitive form, so a centred cell given the way it is published no longer makes the run report that it found no lattice.
* Each reflection is corrected for the sensor's quantum efficiency at the angle it meets the detector (attenuation lengths from the NIST tables, which also fixes the spot-width parallax term on CdTe) and for the attenuation of the flight path between the sample and its pixel; `--flight-path air|helium|vacuum` declares the medium - default air, since no file states it - and the report says what was assumed and what it was worth. The unmerged MTZ records the factors in new `QE` and `FLIGHT` columns beside `LP`, so raw counts are `I / LP * QE * FLIGHT`, and `_process.h5` in new optional `qe` and `flight` datasets.
* Rotation geometry post-refinement fits the crystal and the detector at once, against the observed spot positions and the observed rocking angles together, so the refined distance depends far less on how wrong the file's distance was.
* A coarsely sliced sweep integrates correctly: partials are joined into one rocking event by angle rather than by frame count, so two crossings of the Ewald sphere are no longer summed into one full, and at 0.5 degrees per image or coarser the per-frame geometry refinement accepts a spot whose miss the exposure's own rotation accounts for.
* `rugnux --mode scale` reports the detector tilt and direct beam of the geometry it re-scaled at, instead of zeros that read as a flat detector, and no longer warns that no image was indexed on a run whose lattice came from its input file.
* Every rotation run that determined a space group and merged reports what the mounting cost: `SPINDLE_LOST_UNIQUE_FRACTION=` is the fraction (0-1) of unique reflections the mounting made unmeasurable under the measured point group, also written to the master as `/entry/MX/spindleLostUniqueFraction` and what the mounting warning fires on; `SPINDLE_SYMMETRY_AXIS_ANGLE_DEG=` / `SPINDLE_SYMMETRY_AXIS_ORDER=` describe the mounting in the `--developer` report.
* Stills and grid scans carry a per-image `spindle_blind_fraction` - how much of a rotation sweep's blind cone this orientation would make unrecoverable, 0.5 and above calling for a second orientation - through the CBOR stream, HDF5 (`/entry/MX/spindleBlindFraction`), the plot and scan-result APIs, and the viewer and frontend plots; an absent value means the frame could not be assessed and is not a 0.
* The results report's `REPORT_VERSION` is 7.

Reviewed-on: #77
Co-authored-by: Filip Leonarski <filip.leonarski@psi.ch>
2026-09-09 07:25:13 +02:00

560 lines
16 KiB
C++

// SPDX-FileCopyrightText: 2024 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
// SPDX-License-Identifier: GPL-3.0-only
#include <cmath>
#include "DatasetSettings.h"
#include "Definitions.h"
#include "JFJochException.h"
#include "CheckPath.h"
#define check_max(param, val, max) if ((val) > (max)) throw JFJochException(JFJochExceptionCategory::InputParameterAboveMax, param)
#define check_min(param, val, min) if ((val) < (min)) throw JFJochException(JFJochExceptionCategory::InputParameterBelowMin, param)
#define check_finite(param, val) if (!std::isfinite(val)) throw JFJochException(JFJochExceptionCategory::InputParameterInvalid, param)
DatasetSettings::DatasetSettings() {
photon_energy_keV = WVL_1A_IN_KEV;
detector_distance_mm = 100;
beam_x_pxl = 0.0;
beam_y_pxl = 0.0;
file_prefix = "test";
ntrigger = 1;
images_per_trigger = 1;
compression = CompressionAlgorithm::BSHUF_LZ4;
data_reduction_factor_serialmx = 1.0;
write_nxmx_hdf5_master = true;
spot_finding_enable = true;
poni_rot_1_rad = 0.0f;
poni_rot_2_rad = 0.0f;
poni_rot_3_rad = 0.0f;
max_spot_count = std::min(MAX_SPOT_COUNT, 250);
detect_ice_rings = false;
}
DatasetSettings &DatasetSettings::ImagesPerTrigger(int64_t input) {
check_max("Total number of images", input, 10*1000*1000);
check_min("Total number of images", input, 0);
images_per_trigger = input;
return *this;
}
DatasetSettings &DatasetSettings::NumTriggers(int64_t input) {
check_max("Total number of triggers", input, 10*1000*1000);
check_min("Total number of triggers", input, 1);
ntrigger = input;
return *this;
}
DatasetSettings &DatasetSettings::PhotonEnergy_keV(float input) {
check_finite("Energy (keV)", input);
check_min("Energy (keV)", input, MIN_ENERGY_KEV);
check_max("Energy (keV)", input, MAX_ENERGY_KEV);
photon_energy_keV = input;
return *this;
}
DatasetSettings &DatasetSettings::BeamX_pxl(float input) {
check_finite("Beam center x", input);
beam_x_pxl = input;
return *this;
}
DatasetSettings &DatasetSettings::BeamY_pxl(float input) {
check_finite("Beam center y", input);
beam_y_pxl = input;
return *this;
}
DatasetSettings &DatasetSettings::DetectorDistance_mm(float input) {
check_finite("Detector distance (mm)", input);
check_min("Detector distance (mm)", input, 1);
detector_distance_mm = input;
return *this;
}
DatasetSettings &DatasetSettings::FilePrefix(std::string input) {
CheckPath(input); // multi-user guard: no absolute path, no '..' traversal
return FilePrefixTrusted(std::move(input));
}
DatasetSettings &DatasetSettings::FilePrefixTrusted(std::string input) {
// Offline/local callers (rugnux, the viewer's processing) supply their own output path, so the
// CheckPath guard applied by FilePrefix() is deliberately skipped here - an absolute path is
// allowed. Not reachable from the broker/writer remote-input path.
if ((input.find("_master.h5") == input.length() - 10) && (input.length() > 10))
file_prefix = input.substr(0, input.length() - 10);
else
file_prefix = input;
return *this;
}
DatasetSettings &DatasetSettings::Compression(CompressionAlgorithm input) {
switch (input) {
case CompressionAlgorithm::NO_COMPRESSION:
case CompressionAlgorithm::BSHUF_LZ4:
case CompressionAlgorithm::BSHUF_ZSTD:
case CompressionAlgorithm::BSHUF_ZSTD_RLE:
case CompressionAlgorithm::BSHUF_ZSTD_RLE_HUFF:
compression = input;
break;
default:
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Invalid value for compression enum parameter");
}
return *this;
}
DatasetSettings &DatasetSettings::SetUnitCell(const std::optional<UnitCell> &cell) {
// If one or more cell parameters are zero - assume cell is not provided
// But if all parameters non-zero - throw exception for non-sense values
if (cell && (cell->a != 0) && (cell->b != 0) && (cell->c != 0)
&& (cell->alpha != 0) && (cell->beta != 0) && (cell->gamma != 0)) {
check_min("Unit cell a", cell->a, 0.1);
check_min("Unit cell b", cell->b, 0.1);
check_min("Unit cell c", cell->c, 0.1);
check_min("Angle alpha", cell->alpha, 5.0);
check_min("Angle beta", cell->beta, 5.0);
check_min("Angle gamma", cell->gamma, 5.0);
check_max("Angle alpha", cell->alpha, 355.0);
check_max("Angle beta", cell->beta, 355.0);
check_max("Angle gamma", cell->gamma, 360);
unit_cell = cell;
} else
unit_cell.reset();
return *this;
}
DatasetSettings &DatasetSettings::SpaceGroupNumber(std::optional<int64_t> input) {
if (input) {
check_min("Space group number", input, 1);
check_max("Space group number", input, 230);
}
// A bare number can only name the reference setting.
if (input)
space_group = *gemmi::find_spacegroup_by_number(static_cast<int>(*input));
else
space_group.reset();
return *this;
}
DatasetSettings &DatasetSettings::SetSpaceGroup(const std::optional<gemmi::SpaceGroup> &input) {
space_group = input;
return *this;
}
DatasetSettings &DatasetSettings::SampleName(std::string input) {
sample_name = input;
return *this;
}
DatasetSettings &DatasetSettings::AttenuatorTransmission(const std::optional<float> &input) {
if (input) {
check_finite("Attenuator transmission", input.value());
check_max("Attenuator transmission", input.value(), 1.0);
check_min("Attenuator transmission", input.value(), 0.0);
}
attenuator_transmission = input;
return *this;
}
DatasetSettings &DatasetSettings::TotalFlux(const std::optional<float> &input) {
if (input) {
check_finite("Total flux", input.value());
check_min("Total flux", input.value(), 0.0);
}
total_flux = input;
return *this;
}
DatasetSettings &DatasetSettings::BeamSizeX_um(const std::optional<float> &input) {
if (input) {
check_finite("Beam size X", input.value());
check_min("Beam size X", input.value(), 0.0);
}
beam_size_x_um = input;
return *this;
}
DatasetSettings &DatasetSettings::BeamSizeY_um(const std::optional<float> &input) {
if (input) {
check_finite("Beam size Y", input.value());
check_min("Beam size Y", input.value(), 0.0);
}
beam_size_y_um = input;
return *this;
}
DatasetSettings &DatasetSettings::Goniometer(const std::optional<GoniometerAxis> &input) {
goniometer = input;
return *this;
}
DatasetSettings &DatasetSettings::HeaderAppendix(const nlohmann::json &input) {
header_appendix = input;
return *this;
}
DatasetSettings &DatasetSettings::ImageAppendix(const nlohmann::json &input) {
image_appendix = input;
return *this;
}
std::optional<float> DatasetSettings::GetAttenuatorTransmission() const {
return attenuator_transmission;
}
std::optional<float> DatasetSettings::GetTotalFlux() const {
return total_flux;
}
std::optional<float> DatasetSettings::GetBeamSizeX_um() const {
return beam_size_x_um;
}
std::optional<float> DatasetSettings::GetBeamSizeY_um() const {
return beam_size_y_um;
}
const std::optional<GoniometerAxis> &DatasetSettings::GetGoniometer() const {
return goniometer;
}
std::optional<GoniometerAxis> &DatasetSettings::Goniometer() {
return goniometer;
}
const nlohmann::json& DatasetSettings::GetHeaderAppendix() const {
return header_appendix;
}
const nlohmann::json& DatasetSettings::GetImageAppendix() const {
return image_appendix;
}
std::optional<UnitCell> DatasetSettings::GetUnitCell() const {
return unit_cell;
}
std::optional<int64_t> DatasetSettings::GetSpaceGroupNumber() const {
if (!space_group)
return std::nullopt;
return space_group->number;
}
const std::optional<gemmi::SpaceGroup> &DatasetSettings::GetSpaceGroup() const {
return space_group;
}
std::string DatasetSettings::GetSampleName() const {
return sample_name;
}
float DatasetSettings::GetPhotonEnergy_keV() const {
return photon_energy_keV;
}
float DatasetSettings::GetBeamX_pxl() const {
return beam_x_pxl;
}
float DatasetSettings::GetBeamY_pxl() const {
return beam_y_pxl;
}
float DatasetSettings::GetDetectorDistance_mm() const {
return detector_distance_mm;
}
Coord DatasetSettings::GetScatteringVector() const {
return {0, 0, photon_energy_keV / WVL_1A_IN_KEV};
}
std::string DatasetSettings::GetFilePrefix() const {
return file_prefix;
}
CompressionAlgorithm DatasetSettings::GetCompressionAlgorithm() const {
return compression;
}
int64_t DatasetSettings::GetNumTriggers() const {
return ntrigger;
}
int64_t DatasetSettings::GetImageNumPerTrigger() const {
return images_per_trigger;
}
DatasetSettings &DatasetSettings::ImagesPerFile(const std::optional<int64_t> &input) {
if (input.has_value())
check_min("Images per file", input.value(), 1);
images_per_file = input;
return *this;
}
std::optional<int64_t> DatasetSettings::GetImagesPerFile() const {
return images_per_file;
}
DatasetSettings &DatasetSettings::LossyCompressionSerialMX(float input) {
check_min("Data reduction factor for serial MX", input, 0.0);
check_max("Data reduction factor for serial MX", input, 1.0);
data_reduction_factor_serialmx = input;
return *this;
}
float DatasetSettings::GetLossyCompressionSerialMX() const {
return data_reduction_factor_serialmx;
}
DatasetSettings &DatasetSettings::RunNumber(const std::optional<uint64_t> &input) {
if (input) {
check_min("Run number", input, 0);
check_max("Run number", input, INT64_MAX);
}
run_number = input;
return *this;
}
DatasetSettings & DatasetSettings::RunName(const std::optional<std::string> &input) {
if (input && input.value().empty())
run_name = {};
else
run_name = input;
return *this;
}
DatasetSettings &DatasetSettings::ExperimentGroup(const std::string &input) {
group = input;
return *this;
}
std::optional<uint64_t> DatasetSettings::GetRunNumber() const {
return run_number;
}
std::optional<std::string> DatasetSettings::GetRunName() const {
return run_name;
}
std::string DatasetSettings::GetExperimentGroup() const {
return group;
}
std::optional<std::chrono::nanoseconds> DatasetSettings::GetImageTime() const {
return image_time;
}
DatasetSettings &DatasetSettings::ImageTime(const std::optional<std::chrono::nanoseconds> &input) {
if (input && (input.value().count() == 0))
image_time = {};
else
image_time = input;
return *this;
}
DatasetSettings &DatasetSettings::LossyCompressionPoisson(std::optional<int64_t> input) {
if (!input || (input == 0))
compression_poisson_factor = {};
else {
check_min("Poisson compression factor", input.value(), 1);
check_max("Poisson compression factor", input.value(), 16);
compression_poisson_factor = input;
}
return *this;
}
std::optional<int64_t> DatasetSettings::GetLossyCompressionPoisson() const {
return compression_poisson_factor;
}
DatasetSettings &DatasetSettings::PixelValueLowThreshold(const std::optional<int64_t> &input) {
if (!input || (input == 0))
pixel_value_low_threshold = {};
else {
check_min("Pixel value low threshold", input.value(), 0);
check_max("Pixel value low threshold", input.value(), INT24_MAX - 1);
pixel_value_low_threshold = input;
}
return *this;
}
DatasetSettings &DatasetSettings::PixelValueHighThreshold(const std::optional<int64_t> &input) {
if (!input || (input == 0))
pixel_value_high_threshold = {};
else {
check_min("Pixel value high threshold", input.value(), 1);
check_max("Pixel value high threshold", input.value(), INT32_MAX);
pixel_value_high_threshold = input;
}
return *this;
}
std::optional<int64_t> DatasetSettings::GetPixelValueLowThreshold() const {
return pixel_value_low_threshold;
}
std::optional<int64_t> DatasetSettings::GetPixelValueHighThreshold() const {
return pixel_value_high_threshold;
}
bool DatasetSettings::IsWriteNXmxHDF5Master() const {
return write_nxmx_hdf5_master;
}
DatasetSettings &DatasetSettings::WriteNXmxHDF5Master(bool input) {
write_nxmx_hdf5_master = input;
return *this;
}
std::optional<bool> DatasetSettings::IsSaveCalibration() const {
return save_calibration;
}
DatasetSettings &DatasetSettings::SaveCalibration(std::optional<bool> input) {
save_calibration = input;
return *this;
}
DatasetSettings &DatasetSettings::GridScan(const std::optional<GridScanSettings> &input) {
grid_scan = input;
return *this;
}
std::optional<GridScanSettings> &DatasetSettings::GridScan() {
return grid_scan;
}
const std::optional<GridScanSettings> &DatasetSettings::GetGridScan() const {
return grid_scan;
}
std::optional<float> DatasetSettings::GetPolarizationFactor() const {
return polarization_factor;
}
DatasetSettings &DatasetSettings::BandwidthFWHM(const std::optional<float> &input) {
bandwidth_fwhm = input;
return *this;
}
std::optional<float> DatasetSettings::GetBandwidthFWHM() const {
return bandwidth_fwhm;
}
float DatasetSettings::GetPoniRot3_rad() const {
return poni_rot_3_rad;
}
float DatasetSettings::GetPoniRot2_rad() const {
return poni_rot_2_rad;
}
float DatasetSettings::GetPoniRot1_rad() const {
return poni_rot_1_rad;
}
DatasetSettings &DatasetSettings::PoniRot1_rad(float input) {
check_finite("PONI rotation 1 (radians)", input);
poni_rot_1_rad = input;
return *this;
}
DatasetSettings &DatasetSettings::PoniRot2_rad(float input) {
check_finite("PONI rotation 2 (radians)", input);
poni_rot_2_rad = input;
return *this;
}
DatasetSettings &DatasetSettings::PoniRot3_rad(float input) {
check_finite("PONI rotation 3 (radians)", input);
poni_rot_3_rad = input;
return *this;
}
DatasetSettings &DatasetSettings::PolarizationFactor(const std::optional<float> &input) {
if (input.has_value()) {
check_finite("Polarization factor", input.value());
check_min("Polarization factor", input.value(), -1.0);
check_max("Polarization factor", input.value(), 1.0);
}
polarization_factor = input;
return *this;
}
std::optional<float> DatasetSettings::GetRingCurrent_mA() const {
return ring_current_mA;
}
DatasetSettings &DatasetSettings::RingCurrent_mA(const std::optional<float> &input) {
if (input.has_value()) {
check_min("Ring current (mA)", input, 0.0);
}
ring_current_mA = input;
return *this;
}
std::optional<float> DatasetSettings::GetSampleTemperature_K() const {
return sample_temperature_K;
}
DatasetSettings &DatasetSettings::SampleTemperature_K(const std::optional<float> &input) {
if (input.has_value()) {
check_min("Sample temperature (K)", input, 0.0);
check_max("Sample temperature (K)", input, 1000.0);
}
sample_temperature_K = input;
return *this;
}
DatasetSettings &DatasetSettings::SpotFindingEnable(bool input) {
spot_finding_enable = input;
return *this;
}
bool DatasetSettings::IsSpotFindingEnabled() const {
return spot_finding_enable;
}
DatasetSettings &DatasetSettings::MaxSpotCount(int64_t input) {
check_min("Max spot count", input, MIN_SPOT_COUNT);
check_max("Max spot count", input, MAX_SPOT_COUNT);
max_spot_count = input;
return *this;
}
DatasetSettings & DatasetSettings::DetectIceRings(bool input) {
detect_ice_rings = input;
return *this;
}
bool DatasetSettings::IsDetectIceRings() const {
return detect_ice_rings;
}
DatasetSettings &DatasetSettings::FluorescenceSpectrum(const XrayFluorescenceSpectrum &input) {
fluorescence_spectrum = input;
return *this;
}
const XrayFluorescenceSpectrum & DatasetSettings::GetFluorescenceSpectrum() const {
return fluorescence_spectrum;
}
int64_t DatasetSettings::GetMaxSpotCount() const {
return max_spot_count;
}
std::optional<SmargonPosition> DatasetSettings::GetSmargonPosition() const {
return smargon_position;
}
DatasetSettings &DatasetSettings::Smargon(const std::optional<SmargonPosition> &input) {
smargon_position = input;
return *this;
}