Three reporting changes. None of them touches a decision: the merged .mtz, .hkl, .cif and _image.dat are byte-identical before and after on two test crystals, and the space group and cell are unchanged. REPORT_VERSION moves to 3 and one warning appears where the second change is working. The alternatives line was misleading. "Best space group: C2 or P21 or P2 (indistinguishable from these data)" sat forty lines below the run's single cell, which is C2's. P2 and P21 are primitive: their lattice is a sub-lattice of that one, their cell has half its volume, and their indices are not the ones in the written files. A user acting on "P21" got a group name with no cell, beside a cell belonging to a different group, with nothing saying so. The report now names the alternatives whose centring differs from the chosen group's, gives the volume ratio of the cell each implies, and says that adopting one means reindexing. It stays byte-identical when every alternative shares the chosen centring, which is every other case that occurs - the enantiomorphic and origin-ambiguous pairs, where one cell really does describe them all. It does not print derived cell constants: the primitive cell of a C-centred monoclinic is not P2's conventional cell, and deriving that is a reduction plus a per-group setting choice, which is more than a report line should carry. A centring the data could not test now says so. Prediction runs in P so the search can confirm a centring from the reflections it extinguishes, but when the indexer returns the primitive sub-cell those reflections are never predicted, and an absent count of 0 reads as "predicts no absences" when it means "none was measured". The search abstains correctly and the centring is then taken from the metric, but none of that reached the user. The candidate table gains a three-valued centring column and the report raises a warning where the metric decides it. No decision logic changed - the flag is two counts the search had already made, and nothing reads it back. The twin-law H is now reported for every operator, and the adopted promotion's ratio against its bound. It was previously visible only inside a refusal message, which is why recalibrating that bound recently required rebuilding the statistic from reference data rather than reading it out of runs we already had. The ratio is reported at the Stage A promotion where it is computed, not in the Stage B table, whose rows all share one point group and would have printed the same number on every line. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01CHMmeM1d489zvNFT7ZMN2P
185 lines
8.9 KiB
C++
185 lines
8.9 KiB
C++
// SPDX-FileCopyrightText: 2026 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
|
|
// SPDX-License-Identifier: GPL-3.0-only
|
|
|
|
#include <catch2/catch_all.hpp>
|
|
|
|
#include "../common/DiffractionExperiment.h"
|
|
#include "../image_analysis/scale_merge/Merge.h"
|
|
#include "../reader/JFJochHDF5Reader.h"
|
|
#include "../rugnux/ResultReport.h"
|
|
#include "../writer/FileWriter.h"
|
|
|
|
// Two synthetic ranges, one of each shape the report has to handle.
|
|
namespace {
|
|
SweepQuality TestSweepQuality() {
|
|
SweepQuality out;
|
|
out.measured = true;
|
|
out.sweep_deg = 60.0f;
|
|
out.ranges.push_back(SweepQualityRange{.first_image = 100, .last_image = 149,
|
|
.reason = SweepQualityReason::CrystalOutOfBeam,
|
|
.severity = 0.83f, .rotation_deg = 5.0f,
|
|
.mean_relative_scale = 0.12f, .mean_relative_cc = 0.30f,
|
|
.indexed_fraction = 0.015f});
|
|
out.ranges.push_back(SweepQualityRange{.first_image = 400, .last_image = 499,
|
|
.reason = SweepQualityReason::RadiationDamage,
|
|
.severity = 0.41f, .rotation_deg = 10.0f,
|
|
.mean_relative_scale = 0.55f, .mean_relative_cc = 0.80f,
|
|
.indexed_fraction = 0.62f});
|
|
return out;
|
|
}
|
|
|
|
std::vector<std::string> ReasonVocabulary() {
|
|
std::vector<std::string> out;
|
|
for (int r = 0; r <= static_cast<int>(SweepQualityReason::RadiationDamage); ++r)
|
|
out.emplace_back(SweepQualityReasonCode(static_cast<SweepQualityReason>(r)));
|
|
return out;
|
|
}
|
|
}
|
|
|
|
TEST_CASE("SweepQuality_ReasonVocabulary", "[Diagnostics]") {
|
|
// The codes are an interface - they are written verbatim into <prefix>_report.txt and into HDF5,
|
|
// where the per-image codes index this list from 1. A rename or a reorder silently breaks every
|
|
// consumer, so pin both the spelling and the order here.
|
|
const auto codes = ReasonVocabulary();
|
|
REQUIRE(codes.size() == 5);
|
|
CHECK(codes[0] == "no_diffraction");
|
|
CHECK(codes[1] == "crystal_out_of_beam");
|
|
CHECK(codes[2] == "weak_diffraction");
|
|
CHECK(codes[3] == "loss_of_centring");
|
|
CHECK(codes[4] == "radiation_damage");
|
|
}
|
|
|
|
TEST_CASE("ResultReport_Render", "[Diagnostics]") {
|
|
DiffractionExperiment x(DetJF(1));
|
|
x.ImagesPerTrigger(600);
|
|
|
|
ProcessResult result;
|
|
result.images_processed = 600;
|
|
result.indexing_rate = 0.87f;
|
|
result.consensus_cell = UnitCell{.a = 79.0f, .b = 79.0f, .c = 38.0f,
|
|
.alpha = 90.0f, .beta = 90.0f, .gamma = 90.0f};
|
|
result.space_group_number = 96;
|
|
result.used_beam_x_pxl = 766.62f;
|
|
result.used_beam_y_pxl = 846.87f;
|
|
result.used_distance_mm = 243.53f;
|
|
result.has_merge_statistics = true;
|
|
result.merge_statistics.sweep_quality = TestSweepQuality();
|
|
|
|
const auto text = RenderResultReport("prefix", "in.h5", x, result);
|
|
|
|
// The stable keys a consumer greps for.
|
|
CHECK(text.find("\nREPORT_VERSION= 3\n") != std::string::npos);
|
|
CHECK(text.find("\nOUTPUT_PREFIX= prefix\n") != std::string::npos);
|
|
CHECK(text.find("\nINDEXING_RATE= 0.8700\n") != std::string::npos);
|
|
CHECK(text.find("\nSPACE_GROUP_NUMBER= 96\n") != std::string::npos);
|
|
CHECK(text.find("\nSWEEP_QUALITY_STATUS= COMPUTED\n") != std::string::npos);
|
|
CHECK(text.find("\nSWEEP_QUALITY_COUNT= 2\n") != std::string::npos);
|
|
CHECK(text.find("\nSWEEP_QUALITY_REASONS= no_diffraction crystal_out_of_beam weak_diffraction "
|
|
"loss_of_centring radiation_damage\n") != std::string::npos);
|
|
|
|
// One table row per range, with the reason code verbatim.
|
|
CHECK(text.find(" 100 149 50 5.0 crystal_out_of_beam ")
|
|
!= std::string::npos);
|
|
CHECK(text.find(" 400 499 100 10.0 radiation_damage ")
|
|
!= std::string::npos);
|
|
|
|
// ... and one plain-English WARNING line per range, greppable by the marker alone.
|
|
CHECK(text.find("\nWARNING_COUNT= 2\n") != std::string::npos);
|
|
CHECK(text.find("\nWARNING: Frames 100-149 out of beam (5.0 deg,") != std::string::npos);
|
|
CHECK(text.find("\nWARNING: Frames 400-499 radiation damage (10.0 deg,") != std::string::npos);
|
|
}
|
|
|
|
TEST_CASE("ResultReport_RenderEmpty", "[Diagnostics]") {
|
|
// A clean run and a run that never looked must be distinguishable: both have a count of 0, and
|
|
// only the STATUS key separates them. This is the property a consumer relies on.
|
|
DiffractionExperiment x(DetJF(1));
|
|
|
|
ProcessResult clean;
|
|
clean.has_merge_statistics = true;
|
|
clean.merge_statistics.sweep_quality.measured = true;
|
|
const auto clean_text = RenderResultReport("p", "in.h5", x, clean);
|
|
CHECK(clean_text.find("\nSWEEP_QUALITY_STATUS= COMPUTED\n") != std::string::npos);
|
|
CHECK(clean_text.find("\nSWEEP_QUALITY_COUNT= 0\n") != std::string::npos);
|
|
CHECK(clean_text.find("FIRST_IMAGE LAST_IMAGE") != std::string::npos);
|
|
|
|
ProcessResult not_merged;
|
|
const auto not_merged_text = RenderResultReport("p", "in.h5", x, not_merged);
|
|
CHECK(not_merged_text.find("\nSWEEP_QUALITY_STATUS= NOT_COMPUTED\n") != std::string::npos);
|
|
CHECK(not_merged_text.find("\nMERGE= NOT_PERFORMED\n") != std::string::npos);
|
|
CHECK(not_merged_text.find("FIRST_IMAGE LAST_IMAGE") != std::string::npos);
|
|
}
|
|
|
|
TEST_CASE("ResultReport_RadiationDamage", "[Diagnostics]") {
|
|
// RADIATION_DAMAGE_RELATIVE_B is a number only when there is one. A curve that was measured but
|
|
// that no straight line describes, and a monitor that could not run at all, are different answers,
|
|
// and a consumer has to be able to tell them apart - and both from a measured zero.
|
|
DiffractionExperiment x(DetJF(1));
|
|
ProcessResult result;
|
|
result.has_merge_statistics = true;
|
|
result.radiation_damage_text = "per-batch relative-B";
|
|
result.merge_statistics.radiation_damage_batch_deg = 10.0;
|
|
result.merge_statistics.radiation_damage_b_batch = {0.0f, 4.0f, NAN};
|
|
|
|
result.merge_statistics.radiation_damage_delta_b = 8.25;
|
|
CHECK(RenderResultReport("p", "in.h5", x, result).find("\nRADIATION_DAMAGE_RELATIVE_B= 8.25\n")
|
|
!= std::string::npos);
|
|
|
|
result.merge_statistics.radiation_damage_delta_b = NAN;
|
|
CHECK(RenderResultReport("p", "in.h5", x, result).find("\nRADIATION_DAMAGE_RELATIVE_B= NOT_A_TREND\n")
|
|
!= std::string::npos);
|
|
|
|
result.merge_statistics.radiation_damage_b_batch.clear();
|
|
CHECK(RenderResultReport("p", "in.h5", x, result).find("\nRADIATION_DAMAGE_RELATIVE_B= NOT_MEASURED\n")
|
|
!= std::string::npos);
|
|
}
|
|
|
|
TEST_CASE("SweepQuality_HDF5RoundTrip", "[HDF5][Full][Diagnostics]") {
|
|
// The per-image codes have to survive the writer and come back out of the reader. Without an
|
|
// assertion here the field can ship as all-zeros without anyone noticing.
|
|
DiffractionExperiment x(DetJF(1));
|
|
x.ImagesPerTrigger(6).Compression(CompressionAlgorithm::NO_COMPRESSION)
|
|
.FilePrefix("sweep_quality_roundtrip");
|
|
x.SetFileWriterFormat(FileWriterFormat::NXmxIntegrated).OverwriteExistingFiles(true);
|
|
|
|
// Images 2-3 out of beam, image 5 dead from radiation damage; the rest in no flagged range.
|
|
const std::vector<uint8_t> expected{0, 0,
|
|
static_cast<uint8_t>(SweepQualityReason::CrystalOutOfBeam) + 1,
|
|
static_cast<uint8_t>(SweepQualityReason::CrystalOutOfBeam) + 1,
|
|
0,
|
|
static_cast<uint8_t>(SweepQualityReason::RadiationDamage) + 1};
|
|
|
|
{
|
|
RegisterHDF5Filter();
|
|
StartMessage start_message;
|
|
x.FillMessage(start_message);
|
|
|
|
EndMessage end_message;
|
|
end_message.max_image_number = x.GetImageNum();
|
|
end_message.sweep_quality = expected;
|
|
end_message.sweep_quality_reasons = ReasonVocabulary();
|
|
|
|
FileWriter writer(start_message);
|
|
std::vector<int16_t> image(x.GetPixelsNum(), 42);
|
|
for (int i = 0; i < x.GetImageNum(); i++) {
|
|
DataMessage message{};
|
|
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
|
|
message.number = i;
|
|
REQUIRE_NOTHROW(writer.Write(message));
|
|
}
|
|
writer.WriteHDF5(end_message);
|
|
writer.Finalize();
|
|
}
|
|
|
|
{
|
|
JFJochHDF5Reader reader;
|
|
REQUIRE_NOTHROW(reader.ReadFile("sweep_quality_roundtrip_master.h5"));
|
|
auto dataset = reader.GetDataset();
|
|
CHECK(dataset->sweep_quality == expected);
|
|
// The vocabulary travels with the codes, so a consumer can name them without this source.
|
|
CHECK(dataset->sweep_quality_reasons == ReasonVocabulary());
|
|
}
|
|
|
|
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
|
|
remove("sweep_quality_roundtrip_master.h5");
|
|
}
|