mx/Jungfraujoch
2
1
Fork 0
Code Issues Pull Requests Actions 2 Packages Projects Releases 73 Wiki Activity
Files
a27c4cf26fc86c4e9922b8c66535ef541967b1c8
Jungfraujoch/tests/JFJochReaderTest.cpp
T
leonarski_fandClaude Opus 5 a27c4cf26f reader: take a miniCBF's mounting from the imgCIF axis table its header states 
A miniCBF header states three things about how the instrument is put together
that the reader was assuming instead: which laboratory direction the image's
columns run along, which its rows run along, and which the spindle turns about.
Some beamlines append a CBF template block holding the full imgCIF axis table,
which says all three outright.

Two instruments in the corpus are not what was assumed, in two different ways.
One mounts its detector a quarter turn round, so the image's columns run
vertically. Another turns its spindle about the VERTICAL, with the image mounted
the usual way; its table says so, and its "# Oscillation_axis" line says so a
second way, by naming the image direction the spindle runs along rather than a
vector. Either error leaves the spindle 90 degrees from the image. That is not a
sign, so the run's axis-sign rescue cannot reach it, and no refinement recovers
it: all three affected sweeps indexed nothing usable.

So the table is read. The element axes give the image orientation, matched against
the eight discrete mountings exactly as the NXmx module directions already are -
the match itself moves to DetectorOrientation, so both readers share one
definition rather than two copies. The goniometer axis with no parent gives the
spindle DIRECTION; its sign stays the rescue's business, which is the part a
convention can legitimately differ on. The detector axis with no parent gives the
2theta arm, replacing the assumption that the arm shares the spindle's axis - the
one header stating both states them with the same vector, so this changes no
answer, only what it rests on. imgCIF's frame differs from the internal one by a
half turn about x, a rotation and not a mirror, as writer/HDF5NXmx.cpp already
records from the other side.

Where a header carries no table, a "+SLOW" on the Oscillation_axis line still
says the spindle runs along the image's slow direction. That is the only thing one
of the three affected sets says about it. The axis NAME on that line stays
unusable - the header that carries both says "X.CW" where its own table says Y -
but the direction token is not: where both are present they agree, which is what
makes reading it evidence rather than a guess.

Also: naming a frame with no directory at all now finds its sweep. parent_path()
of a bare filename is empty and iterating an empty path finds nothing, so running
from inside the data directory reported that no images were found.

Measured, with nothing on the command line. The vertical-spindle protein set goes
from no usable lattice to 100% indexed, P 6(3) 2 2 with a cell 0.43% from
deposited, 87846 reflections at 86.3% completeness and CC(1/2) 0.995. Its
companion from the same detector, which has no table and only the +SLOW token,
goes from a spurious monoclinic cell at 2.3% completeness and I/sigma 0.21 to the
right orthorhombic lattice, 97.7% indexed, 59.7% complete, CC(1/2) 0.996. The
quarter-turned set's three sweeps, at three arm positions, now all index without
the hand-passed quarter turn they needed and agree on one cell to 0.03 A. Six
miniCBF sets that state no table and no +SLOW - including one whose
Oscillation_axis line names an axis in a third dialect - are byte-identical in
.hkl, .mtz, .cif and the image statistics, as are two NXmx sets, which is the
shared orientation matcher moving nothing on that path either.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01T3yNBXk4wKdMZy1ak2NY7f
2026-08-30 14:03:32 +02:00

3938 lines
160 KiB
C++
Raw Blame History

// SPDX-FileCopyrightText: 2025 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 "../common/ScanResultGenerator.h"
#include "../writer/FileWriter.h"
#include "../reader/JFJochHDF5Reader.h"
#include "../reader/JFJochCBFReader.h"
#include "../reader/MiniCBF.h"
#include "../compression/JFJochCompressor.h"
#include <fstream>
#include <future>
#include <iomanip>
#include <sstream>
TEST_CASE("HDF5DataType_Sign","[HDF5]") {
HDF5DataType type_u8((uint8_t)0), type_fl(0.0f), type_i32((int32_t) 0), type_u32((uint32_t) 0);
CHECK(!type_u8.IsSigned());
CHECK(type_fl.IsSigned());
CHECK(type_i32.IsSigned());
CHECK(!type_u32.IsSigned());
}
TEST_CASE("HDF5DataType_ElemSize","[HDF5]") {
HDF5DataType type_u8((uint8_t)0), type_fl(0.0f), type_i32((int32_t) 0), type_u32((uint32_t) 0);
CHECK(type_u8.GetElemSize() == 1);
CHECK(type_fl.GetElemSize() == 4);
CHECK(type_i32.GetElemSize() == 4);
CHECK(type_u32.GetElemSize() == 4);
}
TEST_CASE("HDF5DataType_ElemType","[HDF5]") {
HDF5DataType type_u8((uint8_t)0), type_fl(0.0f), type_i32((int32_t) 0), type_u32((uint32_t) 0);
CHECK(type_u8.IsInteger());
CHECK(!type_fl.IsInteger());
CHECK(type_fl.IsFloat());
CHECK(type_i32.IsInteger());
CHECK(type_u32.IsInteger());
}
TEST_CASE("JFJochReader_MasterFile", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test08").ImagesPerTrigger(950).OverwriteExistingFiles(true);
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10))
.SetUnitCell(UnitCell{.a= 10, .b= 20, .c= 30, .alpha= 90, .beta= 101, .gamma = 90});
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
EndMessage end_message;
end_message.max_image_number = 0;
std::unique_ptr<NXmx> master = std::make_unique<NXmx>(start_message);
master->Finalize(end_message);
master.reset();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test08_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetBeamX_pxl() == Catch::Approx(x.GetBeamX_pxl()));
CHECK(dataset->experiment.GetBeamY_pxl() == Catch::Approx(x.GetBeamY_pxl()));
CHECK(dataset->experiment.GetDetectorDistance_mm() == Catch::Approx(x.GetDetectorDistance_mm()));
CHECK(dataset->experiment.GetFrameTime() == x.GetFrameTime());
CHECK(dataset->experiment.GetFrameCountTime() == x.GetFrameCountTime());
CHECK(dataset->experiment.GetWavelength_A() == Catch::Approx(x.GetWavelength_A()));
CHECK(dataset->experiment.GetImageNum() == 0);
REQUIRE(dataset->experiment.GetUnitCell().has_value());
CHECK(dataset->experiment.GetUnitCell()->b == 20.0);
CHECK(dataset->experiment.GetUnitCell()->beta == 101.0);
CHECK(dataset->calibration_data.empty());
}
remove("test08_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_MasterFile_Calibration", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_reader_calibration").ImagesPerTrigger(1).OverwriteExistingFiles(true);
RegisterHDF5Filter();
std::vector<uint16_t> calib_1(200*300, 10);
std::vector<int32_t> calib_2(100*400, 55);
std::vector<float> calib_f(100*400, 1234.56f);
JFJochBitShuffleCompressor compressor(CompressionAlgorithm::BSHUF_LZ4);
auto calib_1_compressed = compressor.Compress(calib_1);
{
StartMessage start_message;
x.FillMessage(start_message);
CompressedImage calibration_01(calib_1, 200, 300);
CompressedImage calibration_02(calib_2, 100, 400);
CompressedImage calibration_f(calib_f, 100, 400);
CompressedImage calibration_01_lz4(
calib_1_compressed.data(), calib_1_compressed.size(),
200, 300, CompressedImageMode::Uint16, CompressionAlgorithm::BSHUF_LZ4
);
calibration_01.Channel("c1");
calibration_02.Channel("c2");
calibration_f.Channel("cf");
calibration_01_lz4.Channel("c1_lz4");
EndMessage end_message;
end_message.max_image_number = 0;
std::unique_ptr<NXmx> master = std::make_unique<NXmx>(start_message);
master->WriteCalibration(calibration_01);
master->WriteCalibration(calibration_01_lz4);
master->WriteCalibration(calibration_02);
master->WriteCalibration(calibration_f);
master->Finalize(end_message);
master.reset();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_reader_calibration_master.h5"));
auto dataset = reader.GetDataset();
REQUIRE(dataset->calibration_data.size() == 4);
CHECK(dataset->calibration_data[0] == "c1");
CHECK(dataset->calibration_data[1] == "c1_lz4");
CHECK(dataset->calibration_data[2] == "c2");
CHECK(dataset->calibration_data[3] == "cf");
std::vector<uint8_t> buffer;
std::vector<uint8_t> buff_2;
REQUIRE_THROWS(reader.ReadCalibration(buffer, "c3"));
CompressedImage test;
REQUIRE_NOTHROW(test = reader.ReadCalibration(buffer, "c1"));
CHECK(test.GetByteDepth() == 2);
CHECK(test.GetHeight() == 300);
CHECK(test.GetWidth() == 200);
CHECK(test.GetMode() == CompressedImageMode::Uint16);
CHECK(reinterpret_cast<const uint16_t *>(test.GetUncompressedPtr(buff_2))[76] == 10);
REQUIRE_NOTHROW(test = reader.ReadCalibration(buffer, "c1_lz4"));
CHECK(test.GetByteDepth() == 2);
CHECK(test.GetHeight() == 300);
CHECK(test.GetWidth() == 200);
CHECK(test.GetMode() == CompressedImageMode::Uint16);
CHECK(reinterpret_cast<const uint16_t *>(test.GetUncompressedPtr(buff_2))[76] == 10);
REQUIRE_NOTHROW(test = reader.ReadCalibration(buffer, "c2"));
CHECK(test.GetByteDepth() == 4);
CHECK(test.GetHeight() == 400);
CHECK(test.GetWidth() == 100);
CHECK(test.GetMode() == CompressedImageMode::Int32);
CHECK(reinterpret_cast<const int32_t *>(test.GetUncompressedPtr(buff_2))[76] == 55);
REQUIRE_NOTHROW(test = reader.ReadCalibration(buffer, "cf"));
CHECK(test.GetByteDepth() == 4);
CHECK(test.GetHeight() == 400);
CHECK(test.GetWidth() == 100);
CHECK(test.GetMode() == CompressedImageMode::Float32);
CHECK(reinterpret_cast<const float *>(test.GetUncompressedPtr(buff_2))[76] == Catch::Approx(1234.56f));
}
remove("test_reader_calibration_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_DefaultExperiment", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_def").OverwriteExistingFiles(true);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
EndMessage end_message;
end_message.max_image_number = 0;
std::unique_ptr<NXmx> master = std::make_unique<NXmx>(start_message);
master->Finalize(end_message);
master.reset();
}
{
JFJochHDF5Reader reader;
DiffractionExperiment x1;
IndexingSettings is;
is.FFT_NumVectors(1024);
x1.ImportIndexingSettings(is);
reader.Experiment(x1);
REQUIRE_NOTHROW(reader.ReadFile("test_def_master.h5"));
auto dataset = reader.GetDataset();
REQUIRE(x1.GetIndexingSettings().GetFFT_NumVectors() == 1024);
}
remove("test_def_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_PixelMask", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test16").ImagesPerTrigger(950).OverwriteExistingFiles(true);
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
RegisterHDF5Filter();
std::vector<uint32_t> pixel_mask(x.GetPixelsNum(), 0);
pixel_mask[5767] = 1;
pixel_mask[x.GetPixelsNum() - 1] = 4;
pixel_mask[0] = 256;
pixel_mask[3] = 1u << PixelMask::BeamStopPixelBit;
ScanResultGenerator generator(x);
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
{
StartMessage start_message;
x.FillMessage(start_message);
start_message.pixel_mask["default"] = pixel_mask;
FileWriter file_set(start_message);
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 0;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
REQUIRE_NOTHROW(generator.Add(message));
EndMessage end_message;
end_message.max_image_number = 1;
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
reader.ReadFile("test16_master.h5");
auto dataset = reader.GetDataset();
REQUIRE(dataset->pixel_mask->GetMask().size() == x.GetPixelsNum());
CHECK(dataset->pixel_mask->GetMask() == pixel_mask);
std::shared_ptr<JFJochReaderImage> reader_image;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(0));
REQUIRE(reader_image);
CHECK(reader_image->Image().at(5767) == GAP_PXL_VALUE);
CHECK(reader_image->Image().at(0) == ERROR_PXL_VALUE);
CHECK(reader_image->Image().at(1) == 0);
CHECK(reader_image->Image().at(2) == 0);
// The beam-stop shadow reads back as its own marker, not as a bad pixel
CHECK(reader_image->Image().at(3) == BEAM_STOP_PXL_VALUE);
CHECK(reader_image->Image().at(x.GetPixelsNum() - 1) == ERROR_PXL_VALUE);
}
remove("test16_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_ROIDefinitions", "[HDF5][Full]") {
RegisterHDF5Filter();
// ROI definitions and the bitmap live in the master file for every format.
auto format = GENERATE(FileWriterFormat::NXmxVDS, FileWriterFormat::NXmxIntegrated);
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_roi").ImagesPerTrigger(950).OverwriteExistingFiles(true)
.SetFileWriterFormat(format);
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
ROIDefinition defs;
defs.boxes.emplace_back("mybox", 10, 20, 30, 40);
defs.circles.emplace_back("mycircle", 100, 200, 15);
defs.azimuthal.emplace_back("mywedge", 2.0f, 4.0f, 30.0f, 90.0f);
x.ROI().SetROI(defs);
ScanResultGenerator generator(x);
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
{
StartMessage start_message;
x.FillMessage(start_message);
start_message.rois = x.ROI().ExportMetadata();
start_message.roi_map = x.ExportROIMap();
FileWriter file_set(start_message);
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 0;
for (const auto &name : {"mybox", "mycircle", "mywedge"})
message.roi[name] = ROIMessage{.sum = 100, .sum_square = 1000, .max_count = 50,
.pixels = 10, .x_weighted = 500, .y_weighted = 600};
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
EndMessage end_message;
end_message.max_image_number = 1;
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_roi_master.h5"));
auto dataset = reader.GetDataset();
const auto &rd = dataset->experiment.ROI().GetROIDefinition();
REQUIRE(rd.boxes.size() == 1);
REQUIRE(rd.circles.size() == 1);
REQUIRE(rd.azimuthal.size() == 1);
CHECK(rd.boxes[0].GetName() == "mybox");
CHECK(rd.boxes[0].GetXMin() == 10);
CHECK(rd.boxes[0].GetXMax() == 20);
CHECK(rd.circles[0].GetName() == "mycircle");
CHECK(rd.circles[0].GetRadius_pxl() == 15.0f);
CHECK(rd.azimuthal[0].GetName() == "mywedge");
CHECK(rd.azimuthal[0].HasPhi());
CHECK(rd.azimuthal[0].GetPhiMin_deg() == 30.0f);
// bitmap read back with the per-pixel footprint and the name->bit index
CHECK(dataset->roi_map.size() == x.GetXPixelsNumConv() * x.GetYPixelsNumConv());
CHECK(dataset->roi_bit_index.size() == 3);
CHECK(dataset->roi_bit_index.at("mybox") == 0);
// per-image ROI results surface from the master (VDS-linked for VDS format)
REQUIRE(dataset->roi.size() == 3);
auto it = std::find(dataset->roi.begin(), dataset->roi.end(), "mybox");
REQUIRE(it != dataset->roi.end());
const size_t idx = std::distance(dataset->roi.begin(), it);
CHECK(dataset->roi_sum.at(idx).at(0) == 100);
}
remove("test_roi_master.h5");
remove("test_roi_data_000001.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_Goniometer", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test17").ImagesPerTrigger(950).OverwriteExistingFiles(true);
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
x.Goniometer(GoniometerAxis("omega", 95, 0.1f, Coord(0,-1,0),{}).ScreeningWedge(0.01f));
RegisterHDF5Filter();
ScanResultGenerator generator(x);
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < 5; i++) {
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
REQUIRE_NOTHROW(generator.Add(message));
}
EndMessage end_message;
end_message.max_image_number = 5;
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
reader.ReadFile("test17_master.h5");
auto dataset = reader.GetDataset();
REQUIRE(!dataset->experiment.GetGridScan().has_value());
REQUIRE(dataset->experiment.GetGoniometer().has_value());
CHECK(dataset->experiment.GetGoniometer()->GetStart_deg() == 95.0);
CHECK(dataset->experiment.GetGoniometer()->GetIncrement_deg() == Catch::Approx(0.1f).margin(0.00001f));
CHECK(dataset->experiment.GetGoniometer()->GetWedge_deg() == Catch::Approx(0.01f).margin(0.00001f));
CHECK(dataset->experiment.GetGoniometer()->GetName() == "omega");
CHECK(dataset->experiment.GetGoniometer()->GetAxis().x == 0);
CHECK(dataset->experiment.GetGoniometer()->GetAxis().y == -1);
CHECK(dataset->experiment.GetGoniometer()->GetAxis().z == 0);
}
remove("test17_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// The axis name is free-form in the API, on the wire and in the writer - tests/CBORTest.cpp round
// trips one literally called "z". The reader used to look only for "omega", so a sweep recorded
// under any other name came back as stills, with nothing to indicate it. This is that case.
TEST_CASE("JFJochReader_Goniometer_NonOmegaName", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test17b").ImagesPerTrigger(950).OverwriteExistingFiles(true);
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
x.Goniometer(GoniometerAxis("phi", 12, 0.2f, Coord(-1,0,0),{}));
RegisterHDF5Filter();
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < 5; i++) {
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = 5;
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
reader.ReadFile("test17b_master.h5");
auto dataset = reader.GetDataset();
REQUIRE(dataset->experiment.GetGoniometer().has_value());
CHECK(dataset->experiment.GetGoniometer()->GetName() == "phi");
CHECK(dataset->experiment.GetGoniometer()->GetStart_deg() == 12.0);
CHECK(dataset->experiment.GetGoniometer()->GetIncrement_deg() == Catch::Approx(0.2f).margin(0.00001f));
CHECK(dataset->experiment.GetGoniometer()->IsScanning());
CHECK(dataset->experiment.GetGoniometer()->GetAxis().x == -1);
}
remove("test17b_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_GridScan", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.PixelSigned(false).BitDepthImage(16).OverwriteExistingFiles(true);
x.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
DatasetSettings d;
d.FilePrefix("test_reader_grid_scan").ImagesPerTrigger(5);
d.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.PhotonEnergy_keV(WVL_1A_IN_KEV)
.GridScan(GridScanSettings(3, -7.5, 8.0, true, true));
x.ImportDatasetSettings(d);
RegisterHDF5Filter();
ScanResultGenerator generator(x);
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < 5; i++) {
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
}
EndMessage end_message;
end_message.max_image_number = 5;
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
reader.ReadFile("test_reader_grid_scan_master.h5");
auto dataset = reader.GetDataset();
// A grid scan carries a stationary spindle: NXmx cannot say "no rotation", and a chain of
// translations alone is not readable (dxtbx raises on it). It must not read back as a sweep.
REQUIRE(dataset->experiment.GetGoniometer().has_value());
CHECK(!dataset->experiment.GetGoniometer()->IsScanning());
CHECK(dataset->experiment.GetGoniometer()->GetIncrement_deg() == 0.0f);
REQUIRE(dataset->experiment.GetGridScan().has_value());
CHECK(dataset->experiment.GetGridScan()->IsSnakeScan());
CHECK(dataset->experiment.GetGridScan()->IsVerticalScan());
CHECK(dataset->experiment.GetGridScan()->GetNFast() == 3);
CHECK(dataset->experiment.GetGridScan()->GetNSlow() == 2);
CHECK(dataset->experiment.GetGridScan()->GetNElem() == 6);
CHECK(dataset->experiment.GetGridScan()->GetGridStepX_um() == Catch::Approx(-7.5));
CHECK(dataset->experiment.GetGridScan()->GetGridStepY_um() == Catch::Approx(8.0));
}
{
// That placeholder spindle must carry one entry per image, not a scalar. Our own reader
// copes with either, so the check has to be on the stored shape: a third-party reader takes
// the image count from the innermost axis of the sample chain, and grid_scan_x/y are
// translations and are passed over - so with a scalar here the whole scan reads as one image.
hid_t file = H5Fopen("test_reader_grid_scan_master.h5", H5F_ACC_RDONLY, H5P_DEFAULT);
REQUIRE(file >= 0);
hid_t omega = H5Dopen2(file, "/entry/sample/transformations/omega", H5P_DEFAULT);
REQUIRE(omega >= 0);
hid_t space = H5Dget_space(omega);
CHECK(H5Sget_simple_extent_ndims(space) == 1);
hsize_t dim = 0;
H5Sget_simple_extent_dims(space, &dim, nullptr);
CHECK(dim == 5);
H5Sclose(space);
H5Dclose(omega);
H5Fclose(file);
}
remove("test_reader_grid_scan_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_DataI16", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test09").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum());
image[0] = INT16_MAX;
image[1] = INT16_MIN;
image[2] = 456;
image[3] = -3456;
ScanResultGenerator generator(x);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
image[5678] = i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.indexing_result = (i % 2 == 0);
message.bkg_estimate = i * 345.6;
message.number = i;
message.profile_radius = 123.09;
generator.Add(message);
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test09_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
REQUIRE(dataset->spot_count.size() == 4);
REQUIRE(dataset->bkg_estimate.size() == 4);
REQUIRE(dataset->profile_radius.size() == 4);
REQUIRE_THROWS(reader.LoadImage(4));
std::shared_ptr<JFJochReaderImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.LoadImage(i));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == SATURATED_PXL_VALUE);
CHECK(reader_image->Image()[1] == ERROR_PXL_VALUE);
CHECK(reader_image->Image()[2] == image[2]);
CHECK(reader_image->Image()[3] == image[3]);
CHECK(reader_image->Image()[5678] == i);
CHECK(dataset->indexing_result[i] == (i % 2 == 0));
CHECK(dataset->bkg_estimate[i] == Catch::Approx(i * 345.6));
CHECK(dataset->profile_radius[i] == Catch::Approx(123.09));
}
}
remove("test09_master.h5");
remove("test09_data_000001.h5");
remove("test09_data_000002.h5");
remove("test09_data_000003.h5");
remove("test09_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_DataI16_OldMasterFormat", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test15").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxLegacy).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum());
image[0] = INT16_MAX;
image[1] = INT16_MIN;
image[2] = 456;
image[3] = -3456;
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
image[5678] = i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.indexing_result = (i % 2 == 0);
message.bkg_estimate = i * 345.6;
message.number = i;
message.profile_radius = 1.64;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test15_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
REQUIRE(dataset->spot_count.size() == 4);
REQUIRE(dataset->bkg_estimate.size() == 4);
REQUIRE(dataset->profile_radius.size() == 4);
REQUIRE_THROWS(reader.LoadImage(4));
std::shared_ptr<JFJochReaderImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.LoadImage(i));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == SATURATED_PXL_VALUE);
CHECK(reader_image->Image()[1] == ERROR_PXL_VALUE);
CHECK(reader_image->Image()[2] == image[2]);
CHECK(reader_image->Image()[3] == image[3]);
CHECK(reader_image->Image()[5678] == i);
CHECK(dataset->profile_radius[i] == Catch::Approx(1.64));
CHECK(dataset->indexing_result[i] == (i % 2 == 0));
CHECK(dataset->bkg_estimate[i] == Catch::Approx(i * 345.6));
}
}
remove("test15_master.h5");
remove("test15_data_000001.h5");
remove("test15_data_000002.h5");
remove("test15_data_000003.h5");
remove("test15_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_DataU16", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test10").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(false)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<uint16_t> image(x.GetPixelsNum());
image[0] = UINT16_MAX;
image[1] = INT16_MAX;
image[2] = 456;
ScanResultGenerator generator(x);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
image[5678] = i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.indexing_result = (i % 2 == 0);
message.bkg_estimate = i * 345.6;
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test10_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
REQUIRE_THROWS(reader.LoadImage(4));
std::shared_ptr<JFJochReaderImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.LoadImage(i));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == SATURATED_PXL_VALUE);
CHECK(reader_image->Image()[1] == INT16_MAX);
CHECK(reader_image->Image()[2] == 456);
CHECK(reader_image->Image()[5678] == i);
CHECK(dataset->indexing_result[i] == (i % 2 == 0));
CHECK(dataset->bkg_estimate[i] == Catch::Approx(i * 345.6));
}
}
remove("test10_master.h5");
remove("test10_data_000001.h5");
remove("test10_data_000002.h5");
remove("test10_data_000003.h5");
remove("test10_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_DataI32", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test11").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(32).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int32_t> image(x.GetPixelsNum());
image[0] = INT32_MAX;
image[1] = INT32_MIN;
image[2] = 456;
ScanResultGenerator generator(x);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
image[5678] = i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test11_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
REQUIRE_THROWS(reader.LoadImage(4));
std::shared_ptr<JFJochReaderImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.LoadImage(i));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == SATURATED_PXL_VALUE);
CHECK(reader_image->Image()[1] == ERROR_PXL_VALUE);
CHECK(reader_image->Image()[2] == 456);
CHECK(reader_image->Image()[5678] == i);
}
}
remove("test11_master.h5");
remove("test11_data_000001.h5");
remove("test11_data_000002.h5");
remove("test11_data_000003.h5");
remove("test11_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_DataU32", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test12").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(32).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(false);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<uint32_t> image(x.GetPixelsNum());
image[0] = UINT32_MAX;
image[1] = static_cast<uint32_t>(INT32_MAX) + 50;
image[2] = 456;
image[3] = INT32_MAX;
image[4] = INT32_MAX - 1;
ScanResultGenerator generator(x);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
image[5678] = i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.number = i;
generator.Add(message);
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test12_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
REQUIRE_THROWS(reader.LoadImage(4));
std::shared_ptr<JFJochReaderImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.LoadImage(i));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == INT32_MAX);
CHECK(reader_image->Image()[1] == INT32_MAX);
CHECK(reader_image->Image()[2] == 456);
CHECK(reader_image->Image()[3] == INT32_MAX);
CHECK(reader_image->Image()[4] == INT32_MAX - 1);
CHECK(reader_image->Image()[5678] == i);
}
}
remove("test12_master.h5");
remove("test12_data_000001.h5");
remove("test12_data_000002.h5");
remove("test12_data_000003.h5");
remove("test12_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_Summation", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test30").ImagesPerTrigger(3).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(3).SetFileWriterFormat(FileWriterFormat::NXmxVDS)
.PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image_1(x.GetPixelsNum(),1);
std::vector<int16_t> image_2(x.GetPixelsNum(),2);
std::vector<int16_t> image_3(x.GetPixelsNum(),3);
image_3[0] = INT16_MAX;
image_2[1] = INT16_MIN;
ScanResultGenerator generator(x);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
std::vector<SpotToSave> spots;
DataMessage message{};
message.spots = spots;
message.image = CompressedImage(image_1, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 0;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
message.image = CompressedImage(image_2, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 1;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
message.image = CompressedImage(image_3, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 2;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test30_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 3);
std::shared_ptr<JFJochReaderImage> reader_image;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(0, 3));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == SATURATED_PXL_VALUE);
CHECK(reader_image->Image()[1] == ERROR_PXL_VALUE);
CHECK(reader_image->Image()[2] == 1 + 2 +3);
CHECK(reader_image->Image()[5678] == 1 + 2 +3);
CHECK(reader_image->Image()[x.GetPixelsNum() - 1] == 1 + 2 +3);
REQUIRE_THROWS(reader.LoadImage(1, 3));
}
remove("test30_master.h5");
remove("test30_data_000001.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_Summation_5", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test31").ImagesPerTrigger(5).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(5).SetFileWriterFormat(FileWriterFormat::NXmxVDS)
.PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image_1(x.GetPixelsNum(),1);
std::vector<int16_t> image_2(x.GetPixelsNum(),2);
std::vector<int16_t> image_3(x.GetPixelsNum(),3);
std::vector<int16_t> image_4(x.GetPixelsNum(),4);
std::vector<int16_t> image_5(x.GetPixelsNum(),5);
image_3[0] = INT16_MAX;
image_2[1] = INT16_MIN;
ScanResultGenerator generator(x);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
std::vector<SpotToSave> spots;
DataMessage message{};
message.spots = spots;
message.image = CompressedImage(image_1, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 0;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
message.image = CompressedImage(image_2, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 1;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
message.image = CompressedImage(image_3, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 2;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
message.image = CompressedImage(image_4, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 3;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
message.image = CompressedImage(image_5, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = 4;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test31_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 5);
std::shared_ptr<JFJochReaderImage> reader_image;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(0, 5));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == SATURATED_PXL_VALUE);
CHECK(reader_image->Image()[1] == ERROR_PXL_VALUE);
CHECK(reader_image->Image()[2] == 1 + 2 + 3 + 4 + 5);
CHECK(reader_image->Image()[5678] == 1 + 2 + 3 + 4 + 5);
CHECK(reader_image->Image()[x.GetPixelsNum() - 1] == 1 + 2 + 3 + 4 + 5);
REQUIRE_THROWS(reader.LoadImage(1, 6));
}
remove("test31_master.h5");
remove("test31_data_000001.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_Azint", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test27").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(false);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
AzimuthalIntegrationSettings azint_settings;
azint_settings.AzimuthalBinCount(4);
x.ImportAzimuthalIntegrationSettings(azint_settings);
// The high-q limit is unset, i.e. "as far as the detector reaches", so read the settings back from
// the experiment, where that has been resolved against the geometry - that is what the bins are.
azint_settings = x.GetAzimuthalIntegrationSettings();
std::vector<uint16_t> image(x.GetPixelsNum());
AzimuthalIntegrationMapping azint(x, PixelMask(x));
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
start_message.az_int_bin_to_q = azint.GetBinToQ();
start_message.az_int_bin_to_phi = azint.GetBinToPhi();
start_message.az_int_q_bin_count = azint.GetQBinCount();
start_message.az_int_phi_bin_count = azint.GetAzimuthalBinCount();
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.number = i;
message.az_int_profile = std::vector<float>(azint_settings.GetBinCount(), 57);
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test27_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
std::shared_ptr<JFJochReaderImage> reader_image;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(0));
REQUIRE(reader_image);
CHECK(reader_image->Dataset().az_int_bin_to_q.size() == azint_settings.GetBinCount());
CHECK(reader_image->Dataset().azimuthal_bins == azint_settings.GetAzimuthalBinCount());
CHECK(reader_image->Dataset().q_bins == azint_settings.GetQBinCount());
REQUIRE(reader_image->ImageData().az_int_profile.size() == azint_settings.GetBinCount());
CHECK(reader_image->ImageData().az_int_profile[23] == 57);
CHECK(reader_image->GetAzInt1D_BinToQ().size() == azint_settings.GetQBinCount());
REQUIRE(reader_image->GetAzInt1D().size() == azint_settings.GetQBinCount());
CHECK(reader_image->GetAzInt1D()[23] == 4 * 57);
}
remove("test27_master.h5");
remove("test27_data_000001.h5");
remove("test27_data_000002.h5");
remove("test27_data_000003.h5");
remove("test27_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_NiggliClass", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test95").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxLegacy).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum());
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
LatticeMessage lm{
.centering = 'F',
.niggli_class = 1,
.crystal_system = gemmi::CrystalSystem::Cubic,
};
DataMessage message{};
message.number = 0;
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.indexing_result = true;
message.indexing_lattice = CrystalLattice(40, 50, 60, 90, 90, 90);
message.lattice_type = lm;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
message.number = 1;
message.indexing_result = false;
message.indexing_lattice = std::nullopt;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
EndMessage end_message;
end_message.max_image_number = 2;
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test95_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 2);
std::shared_ptr<JFJochReaderImage> reader_image, reader_image_2;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(0));
REQUIRE(reader_image);
CHECK(reader_image->ImageData().indexing_result.value() == true);
REQUIRE(reader_image->ImageData().indexing_lattice);
REQUIRE(reader_image->ImageData().lattice_type);
CHECK(reader_image->ImageData().lattice_type->centering == 'F');
CHECK(reader_image->ImageData().lattice_type->niggli_class == 1);
CHECK(reader_image->ImageData().lattice_type->crystal_system == gemmi::CrystalSystem::Cubic);
REQUIRE_NOTHROW(reader_image_2 = reader.LoadImage(1));
REQUIRE(reader_image_2);
CHECK(!reader_image_2->ImageData().indexing_result.value());
REQUIRE(!reader_image_2->ImageData().indexing_lattice);
REQUIRE(!reader_image_2->ImageData().lattice_type);
}
remove("test95_master.h5");
remove("test95_data_000001.h5");
remove("test95_data_000002.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_NiggliClass_VDS", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test95").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum());
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
LatticeMessage lm{
.centering = 'F',
.niggli_class = 1,
.crystal_system = gemmi::CrystalSystem::Cubic,
};
DataMessage message_0{};
message_0.number = 0;
message_0.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message_0.indexing_result = true;
message_0.indexing_lattice = CrystalLattice(40, 50, 60, 90, 90, 90);
message_0.lattice_type = lm;
REQUIRE_NOTHROW(file_set.WriteHDF5(message_0));
DataMessage message_1{};
message_1.number = 1;
message_1.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message_1.indexing_result = false;
message_1.indexing_lattice = std::nullopt;
REQUIRE_NOTHROW(file_set.WriteHDF5(message_1));
EndMessage end_message;
end_message.max_image_number = 2;
end_message.image_indexed = {true, false};
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test95_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 2);
std::shared_ptr<JFJochReaderImage> reader_image, reader_image_2;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(0));
REQUIRE(reader_image);
REQUIRE(reader_image->ImageData().indexing_result.has_value());
CHECK(reader_image->ImageData().indexing_result.value() == true);
REQUIRE(reader_image->ImageData().indexing_lattice.has_value());
REQUIRE(reader_image->ImageData().lattice_type.has_value());
CHECK(reader_image->ImageData().lattice_type->centering == 'F');
CHECK(reader_image->ImageData().lattice_type->niggli_class == 1);
CHECK(reader_image->ImageData().lattice_type->crystal_system == gemmi::CrystalSystem::Cubic);
REQUIRE_NOTHROW(reader_image_2 = reader.LoadImage(1));
REQUIRE(reader_image_2);
CHECK(!reader_image_2->ImageData().indexing_result.value());
REQUIRE(!reader_image_2->ImageData().indexing_lattice);
REQUIRE(!reader_image_2->ImageData().lattice_type);
}
remove("test95_master.h5");
remove("test95_data_000001.h5");
remove("test95_data_000002.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_MissingEntries", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test96").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(10).SetFileWriterFormat(FileWriterFormat::NXmxLegacy).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum());
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
message.number = 0;
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.indexing_result = true;
message.indexing_lattice = CrystalLattice(40, 50, 60, 90, 90, 90);
message.spot_count_indexed = 56;
message.spot_count = 85;
message.b_factor = 123.45;
message.spots = {SpotToSave{.x = 10, .y=50, .intensity = 80}};
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
message.number = 1;
message.indexing_result = false;
message.indexing_lattice = std::nullopt;
message.spot_count_indexed = std::nullopt;
message.spot_count = 70;
message.b_factor = std::nullopt;
message.spots = {SpotToSave{.x = 10, .y=50, .intensity = 80}};
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
EndMessage end_message;
end_message.max_image_number = 2;
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test96_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 2);
REQUIRE(dataset->b_factor.size() == 2);
REQUIRE(dataset->spot_count_indexed.size() == 2);
CHECK(dataset->b_factor[0] == Catch::Approx(123.45));
CHECK(std::isnan(dataset->b_factor[1]));
CHECK(dataset->spot_count_indexed[0] == 56);
CHECK(dataset->spot_count_indexed[1] == 0);
std::shared_ptr<JFJochReaderImage> reader_image, reader_image_2;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(0));
REQUIRE(reader_image);
REQUIRE(reader_image->ImageData().b_factor.has_value());
CHECK(reader_image->ImageData().b_factor.value() == Catch::Approx(123.45));
REQUIRE(reader_image->ImageData().spot_count_indexed.has_value());
CHECK(reader_image->ImageData().spot_count_indexed.value() == 56);
REQUIRE_NOTHROW(reader_image_2 = reader.LoadImage(1));
REQUIRE(reader_image_2);
CHECK(reader_image_2->ImageData().spot_count_indexed.has_value());
CHECK(reader_image_2->ImageData().spot_count_indexed.value() == 0);
}
remove("test96_master.h5");
remove("test96_data_000001.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_Spots_OldMasterFormat", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test40").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxLegacy).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum());
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
spots.push_back(SpotToSave{
.x = 1, .y = 2, .intensity = 376,
.h = 11, .k = -3, .l = -5,
.dist_ewald_sphere = 0.1234f,
.ice_ring = true,
.indexed = true
});
spots.push_back(SpotToSave{
.x = 7, .y = -3, .intensity = 0.156f,
.ice_ring = false,
.indexed = false,
});
image[5678] = i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.indexing_result = (i % 2 == 0);
message.number = i;
message.spot_count = 72 + i;
message.spot_count_ice_rings = 45 + 2 * i;
message.spot_count_low_res = 12 + 3 * i;
message.spot_count_indexed = 15 + 4 * i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test40_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
CHECK(dataset->spot_count[1] == 72 + 1);
CHECK(dataset->spot_count_ice_rings[2] == 45 + 2 * 2);
CHECK(dataset->spot_count_low_res[3] == 12 + 3 * 3);
CHECK(dataset->spot_count_indexed[0] == 15);
REQUIRE_THROWS(reader.LoadImage(4));
std::shared_ptr<JFJochReaderImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.LoadImage(i));
REQUIRE(reader_image);
CHECK(reader_image->ImageData().spot_count == 72 + i);
CHECK(reader_image->ImageData().spot_count_ice_rings == 45 + 2 * i);
CHECK(reader_image->ImageData().spot_count_low_res == 12 + 3 * i);
CHECK(reader_image->ImageData().spot_count_indexed == 15 + 4 * i);
REQUIRE(reader_image->ImageData().spots.size() == 2);
CHECK(reader_image->ImageData().spots[0].x == 1);
CHECK(reader_image->ImageData().spots[0].y == 2);
CHECK(reader_image->ImageData().spots[0].intensity == 376);
CHECK(reader_image->ImageData().spots[0].ice_ring == true);
CHECK(reader_image->ImageData().spots[0].indexed == true);
CHECK(reader_image->ImageData().spots[0].h == 11);
CHECK(reader_image->ImageData().spots[0].k == -3);
CHECK(reader_image->ImageData().spots[0].l == -5);
CHECK(reader_image->ImageData().spots[0].dist_ewald_sphere == Catch::Approx(0.1234f));
CHECK(reader_image->ImageData().spots[1].x == 7);
CHECK(reader_image->ImageData().spots[1].y == -3);
CHECK(reader_image->ImageData().spots[1].intensity == Catch::Approx(0.156f));
CHECK(reader_image->ImageData().spots[1].ice_ring == false);
CHECK(reader_image->ImageData().spots[1].indexed == false);
}
}
remove("test40_master.h5");
remove("test40_data_000001.h5");
remove("test40_data_000002.h5");
remove("test40_data_000003.h5");
remove("test40_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_Spots_VDS", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test41").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxLegacy).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum());
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
std::vector<SpotToSave> spots;
spots.push_back(SpotToSave{
.x = 1, .y = 2, .intensity = 376,
.h = 11, .k = -3, .l = -5,
.dist_ewald_sphere = 0.1234f,
.ice_ring = true,
.indexed = true
});
spots.push_back(SpotToSave{
.x = 7, .y = -3, .intensity = 0.156f,
.ice_ring = false,
.indexed = false,
});
image[5678] = i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spots = spots;
message.indexing_result = (i % 2 == 0);
message.number = i;
message.spot_count = 72 + i;
message.spot_count_ice_rings = 45 + 2 * i;
message.spot_count_low_res = 12 + 3 * i;
message.spot_count_indexed = 15 + 4 * i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test41_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
CHECK(dataset->spot_count[1] == 72 + 1);
CHECK(dataset->spot_count_ice_rings[2] == 45 + 2 * 2);
CHECK(dataset->spot_count_low_res[3] == 12 + 3 * 3);
CHECK(dataset->spot_count_indexed[0] == 15);
REQUIRE_THROWS(reader.LoadImage(4));
std::shared_ptr<JFJochReaderImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.LoadImage(i));
REQUIRE(reader_image);
CHECK(reader_image->ImageData().spot_count == 72 + i);
CHECK(reader_image->ImageData().spot_count_ice_rings == 45 + 2 * i);
CHECK(reader_image->ImageData().spot_count_low_res == 12 + 3 * i);
CHECK(reader_image->ImageData().spot_count_indexed == 15 + 4 * i);
REQUIRE(reader_image->ImageData().spots.size() == 2);
CHECK(reader_image->ImageData().spots[0].x == 1);
CHECK(reader_image->ImageData().spots[0].y == 2);
CHECK(reader_image->ImageData().spots[0].intensity == 376);
CHECK(reader_image->ImageData().spots[0].ice_ring == true);
CHECK(reader_image->ImageData().spots[0].indexed == true);
CHECK(reader_image->ImageData().spots[0].h == 11);
CHECK(reader_image->ImageData().spots[0].k == -3);
CHECK(reader_image->ImageData().spots[0].l == -5);
CHECK(reader_image->ImageData().spots[0].dist_ewald_sphere == Catch::Approx(0.1234f));
CHECK(reader_image->ImageData().spots[1].x == 7);
CHECK(reader_image->ImageData().spots[1].y == -3);
CHECK(reader_image->ImageData().spots[1].intensity == Catch::Approx(0.156f));
CHECK(reader_image->ImageData().spots[1].ice_ring == false);
CHECK(reader_image->ImageData().spots[1].indexed == false);
}
}
remove("test41_master.h5");
remove("test41_data_000001.h5");
remove("test41_data_000002.h5");
remove("test41_data_000003.h5");
remove("test41_data_000004.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_InstrumentMetadata_Sample_RingCurrent", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
// Set identifying names and ring current (mA in API; writer stores A)
x.FilePrefix("test_meta").ImagesPerTrigger(0).OverwriteExistingFiles(true);
InstrumentMetadata metadata;
metadata.InstrumentName("PXI").SourceName("SLS");
x.ImportInstrumentMetadata(metadata).SampleName("test_sample").RingCurrent_mA(399.5); // 0.3995 A
x.TotalFlux(1e7).AttenuatorTransmission(0.56);
x.DetectIceRings(false);
// Minimal other required fields
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
// also set fluorescence spectrum
x.FluorescenceSpectrum(XrayFluorescenceSpectrum({1.0f, 2.0f, 3.0f}, {5.0f, 7.0f, 6.0f}));
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
EndMessage end_message;
end_message.max_image_number = 0;
std::unique_ptr<NXmx> master = std::make_unique<NXmx>(start_message);
master->Finalize(end_message);
master.reset();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_meta_master.h5"));
auto dataset = reader.GetDataset();
auto meta = dataset->experiment.GetInstrumentMetadata();
CHECK(meta.GetInstrumentName() == "PXI");
CHECK(meta.GetSourceName() == "SLS");
CHECK(dataset->experiment.GetAttenuatorTransmission() == Catch::Approx(0.56));
CHECK(dataset->experiment.GetTotalFlux() == Catch::Approx(1e7));
// Sample name
CHECK(dataset->experiment.GetSampleName() == "test_sample");
// Ring current read back in mA; allow small fp tolerance
CHECK(dataset->experiment.GetRingCurrent_mA().has_value());
CHECK(dataset->experiment.GetRingCurrent_mA().value() == Catch::Approx(399.5));
CHECK(!dataset->experiment.IsDetectIceRings());
// Fluorescence spectrum presence and values
REQUIRE(!dataset->experiment.GetFluorescenceSpectrum().empty());
const auto& fl = dataset->experiment.GetFluorescenceSpectrum();
CHECK(fl.GetEnergy_eV().size() == 3);
CHECK(fl.GetData().size() == 3);
CHECK(fl.GetEnergy_eV()[0] == Catch::Approx(1.0f));
CHECK(fl.GetData()[1] == Catch::Approx(7.0f));
}
remove("test_meta_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_NXmxIntegrated", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_reader_integrated").ImagesPerTrigger(3).OverwriteExistingFiles(true);
x.BitDepthImage(16).PixelSigned(false).SetFileWriterFormat(FileWriterFormat::NXmxIntegrated);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
AzimuthalIntegrationSettings azint_settings;
azint_settings.AzimuthalBinCount(4);
x.ImportAzimuthalIntegrationSettings(azint_settings);
// The high-q limit is unset, i.e. "as far as the detector reaches", so read the settings back from
// the experiment, where that has been resolved against the geometry - that is what the bins are.
azint_settings = x.GetAzimuthalIntegrationSettings();
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
image[0] = UINT16_MAX;
image[1] = 123;
image[5678] = 321;
AzimuthalIntegrationMapping azint(x, PixelMask(x));
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
start_message.az_int_bin_to_q = azint.GetBinToQ();
start_message.az_int_bin_to_phi = azint.GetBinToPhi();
start_message.az_int_q_bin_count = azint.GetQBinCount();
start_message.az_int_phi_bin_count = azint.GetAzimuthalBinCount();
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
image[5678] = 321 + i;
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
message.image_collection_efficiency = 0.9f + 0.01f * i;
message.az_int_profile = std::vector<float>(azint_settings.GetBinCount(), static_cast<float>(50 + i));
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_reader_integrated_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 3);
REQUIRE(dataset->efficiency.size() == 3);
CHECK(dataset->efficiency[0] == Catch::Approx(0.90f));
CHECK(dataset->efficiency[1] == Catch::Approx(0.91f));
CHECK(dataset->efficiency[2] == Catch::Approx(0.92f));
CHECK(dataset->az_int_bin_to_q.size() == azint_settings.GetBinCount());
CHECK(dataset->azimuthal_bins == azint_settings.GetAzimuthalBinCount());
CHECK(dataset->q_bins == azint_settings.GetQBinCount());
std::shared_ptr<JFJochReaderImage> reader_image;
REQUIRE_NOTHROW(reader_image = reader.LoadImage(1));
REQUIRE(reader_image);
CHECK(reader_image->Image()[0] == SATURATED_PXL_VALUE);
CHECK(reader_image->Image()[1] == 123);
CHECK(reader_image->Image()[5678] == 322);
REQUIRE(reader_image->ImageData().image_collection_efficiency.has_value());
CHECK(reader_image->ImageData().image_collection_efficiency.value() == Catch::Approx(0.91f));
REQUIRE(reader_image->ImageData().az_int_profile.size() == azint_settings.GetBinCount());
CHECK(reader_image->ImageData().az_int_profile[0] == Catch::Approx(51.0f));
CHECK(reader_image->ImageData().az_int_profile[23] == Catch::Approx(51.0f));
}
remove("test_reader_integrated_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_GetRawImage_NXmxLegacy", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_read_raw_image").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(2).SetFileWriterFormat(FileWriterFormat::NXmxLegacy).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::BSHUF_ZSTD);
std::vector<int16_t> image(x.GetPixelsNum());
for (int i = 0; i < image.size(); i++)
image[i] = static_cast<int16_t>((i * 7 + 33) % UINT16_MAX);
RegisterHDF5Filter();
JFJochBitShuffleCompressor compressor(CompressionAlgorithm::BSHUF_ZSTD);
auto compressed_image = compressor.Compress(image);
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(compressed_image, x.GetXPixelsNum(), x.GetYPixelsNum(),
CompressedImageMode::Int16, CompressionAlgorithm::BSHUF_ZSTD);
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_read_raw_image_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
std::shared_ptr<JFJochReaderRawImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.GetRawImage(i));
CHECK(reader_image->image.GetMode() == CompressedImageMode::Int16);
CHECK(reader_image->image.GetCompressionAlgorithm() == CompressionAlgorithm::BSHUF_ZSTD);
CHECK(reader_image->image.GetWidth() == x.GetXPixelsNum());
CHECK(reader_image->image.GetHeight() == x.GetYPixelsNum());
CHECK(reader_image->image.GetCompressedSize() == compressed_image.size());
CHECK(reader_image->image.GetCompressed() == reader_image->image_buffer.data());
REQUIRE(reader_image->image_buffer.size() == compressed_image.size());
CHECK(memcmp(reader_image->image_buffer.data(), compressed_image.data(), compressed_image.size()) == 0);
}
}
remove("test_read_raw_image_master.h5");
remove("test_read_raw_image_data_000001.h5");
remove("test_read_raw_image_data_000002.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// GetRawImage takes the chunk address under the HDF5 lock, then reads the bytes outside it, so this
// is the one path where several workers are inside the reader at once - which is how rugnux uses it.
// The per-image cases above are all single-threaded and would not notice the file being pulled from
// under a read, nor a cache entry racing its own creation.
TEST_CASE("JFJochReader_GetRawImage_Concurrent", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_raw_concurrent").ImagesPerTrigger(16).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(4).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true);
x.Compression(CompressionAlgorithm::BSHUF_ZSTD);
std::vector<int16_t> image(x.GetPixelsNum());
for (size_t i = 0; i < image.size(); i++)
image[i] = static_cast<int16_t>((i * 11 + 5) % UINT16_MAX);
RegisterHDF5Filter();
JFJochBitShuffleCompressor compressor(CompressionAlgorithm::BSHUF_ZSTD);
const auto compressed_image = compressor.Compress(image);
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(compressed_image, x.GetXPixelsNum(), x.GetYPixelsNum(),
CompressedImageMode::Int16, CompressionAlgorithm::BSHUF_ZSTD);
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_raw_concurrent_master.h5"));
std::vector<std::future<bool>> workers;
workers.reserve(8);
for (int w = 0; w < 8; w++) {
workers.push_back(std::async(std::launch::async, [&reader, &compressed_image, &x]() {
for (int i = 0; i < x.GetImageNum(); i++) {
auto raw = reader.GetRawImage(i);
if (raw->image_buffer.size() != compressed_image.size())
return false;
if (memcmp(raw->image_buffer.data(), compressed_image.data(),
compressed_image.size()) != 0)
return false;
}
return true;
}));
}
for (auto &worker: workers)
CHECK(worker.get());
}
remove("test_raw_concurrent_master.h5");
for (int f = 1; f <= 4; f++)
remove(("test_raw_concurrent_data_00000" + std::to_string(f) + ".h5").c_str());
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_GetRawImage_VDS", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_read_raw_image").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(2).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::BSHUF_ZSTD);
std::vector<int16_t> image(x.GetPixelsNum());
for (int i = 0; i < image.size(); i++)
image[i] = static_cast<int16_t>((i * 7 + 33) % UINT16_MAX);
RegisterHDF5Filter();
JFJochBitShuffleCompressor compressor(CompressionAlgorithm::BSHUF_ZSTD);
auto compressed_image = compressor.Compress(image);
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(compressed_image, x.GetXPixelsNum(), x.GetYPixelsNum(),
CompressedImageMode::Int16, CompressionAlgorithm::BSHUF_ZSTD);
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_read_raw_image_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
std::shared_ptr<JFJochReaderRawImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.GetRawImage(i));
CHECK(reader_image->image.GetMode() == CompressedImageMode::Int16);
CHECK(reader_image->image.GetCompressionAlgorithm() == CompressionAlgorithm::BSHUF_ZSTD);
CHECK(reader_image->image.GetWidth() == x.GetXPixelsNum());
CHECK(reader_image->image.GetHeight() == x.GetYPixelsNum());
CHECK(reader_image->image.GetCompressedSize() == compressed_image.size());
CHECK(reader_image->image.GetCompressed() == reader_image->image_buffer.data());
REQUIRE(reader_image->image_buffer.size() == compressed_image.size());
CHECK(memcmp(reader_image->image_buffer.data(), compressed_image.data(), compressed_image.size()) == 0);
}
}
remove("test_read_raw_image_master.h5");
remove("test_read_raw_image_data_000001.h5");
remove("test_read_raw_image_data_000002.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_GetRawImage_Integrated", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_read_raw_image").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(2).SetFileWriterFormat(FileWriterFormat::NXmxIntegrated).PixelSigned(true)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT);
x.Compression(CompressionAlgorithm::BSHUF_ZSTD);
std::vector<int16_t> image(x.GetPixelsNum());
for (int i = 0; i < image.size(); i++)
image[i] = static_cast<int16_t>((i * 7 + 33) % UINT16_MAX);
RegisterHDF5Filter();
JFJochBitShuffleCompressor compressor(CompressionAlgorithm::BSHUF_ZSTD);
auto compressed_image = compressor.Compress(image);
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(compressed_image, x.GetXPixelsNum(), x.GetYPixelsNum(),
CompressedImageMode::Int16, CompressionAlgorithm::BSHUF_ZSTD);
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_read_raw_image_master.h5"));
auto dataset = reader.GetDataset();
CHECK(dataset->experiment.GetImageNum() == 4);
std::shared_ptr<JFJochReaderRawImage> reader_image;
for (int i = 0; i < 4; i++) {
REQUIRE_NOTHROW(reader_image = reader.GetRawImage(i));
CHECK(reader_image->image.GetMode() == CompressedImageMode::Int16);
CHECK(reader_image->image.GetCompressionAlgorithm() == CompressionAlgorithm::BSHUF_ZSTD);
CHECK(reader_image->image.GetWidth() == x.GetXPixelsNum());
CHECK(reader_image->image.GetHeight() == x.GetYPixelsNum());
CHECK(reader_image->image.GetCompressedSize() == compressed_image.size());
CHECK(reader_image->image.GetCompressed() == reader_image->image_buffer.data());
REQUIRE(reader_image->image_buffer.size() == compressed_image.size());
CHECK(memcmp(reader_image->image_buffer.data(), compressed_image.data(), compressed_image.size()) == 0);
}
}
remove("test_read_raw_image_master.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_HDF5DataSource_Integrated", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("source_integrated").ImagesPerTrigger(5).OverwriteExistingFiles(true);
x.BitDepthImage(16).SetFileWriterFormat(FileWriterFormat::NXmxIntegrated).PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 17);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter writer(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
image[5678] = static_cast<int16_t>(100 + i);
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(writer.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
writer.WriteHDF5(end_message);
writer.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("source_integrated_master.h5"));
auto source = reader.GetHDF5DataSource(1, 3);
REQUIRE(source.size() == 1);
CHECK(source[0].filename == "source_integrated_master.h5");
CHECK(source[0].dataset == "/entry/data/data");
CHECK(source[0].source_first_image == 1);
CHECK(source[0].virtual_first_image == 0);
CHECK(source[0].image_count == 3);
}
remove("source_integrated_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_HDF5DataSource_VDS", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("source_vds_mapping").ImagesPerTrigger(5).ImagesPerFile(2).OverwriteExistingFiles(true);
x.BitDepthImage(16).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 21);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter writer(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
image[5678] = static_cast<int16_t>(200 + i);
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(writer.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
writer.WriteHDF5(end_message);
writer.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("source_vds_mapping_master.h5"));
// Range crosses file boundary:
// global images 1,2,3 map to:
// data_000001 image 1
// data_000002 images 0,1
auto source = reader.GetHDF5DataSource(1, 3);
REQUIRE(source.size() == 2);
CHECK(source[0].filename == "source_vds_mapping_data_000001.h5");
CHECK(source[0].dataset == "/entry/data/data");
CHECK(source[0].source_first_image == 1);
CHECK(source[0].virtual_first_image == 0);
CHECK(source[0].image_count == 1);
CHECK(source[1].filename == "source_vds_mapping_data_000002.h5");
CHECK(source[1].dataset == "/entry/data/data");
CHECK(source[1].source_first_image == 0);
CHECK(source[1].virtual_first_image == 1);
CHECK(source[1].image_count == 2);
}
remove("source_vds_mapping_master.h5");
remove("source_vds_mapping_data_000001.h5");
remove("source_vds_mapping_data_000002.h5");
remove("source_vds_mapping_data_000003.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_HDF5DataSource_Legacy", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("source_legacy_mapping").ImagesPerTrigger(5).ImagesPerFile(2).OverwriteExistingFiles(true);
x.BitDepthImage(16).SetFileWriterFormat(FileWriterFormat::NXmxLegacy).PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 31);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter writer(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
image[5678] = static_cast<int16_t>(300 + i);
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(writer.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
writer.WriteHDF5(end_message);
writer.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("source_legacy_mapping_master.h5"));
auto source = reader.GetHDF5DataSource(1, 3);
REQUIRE(source.size() == 2);
CHECK(source[0].filename == "source_legacy_mapping_data_000001.h5");
CHECK(source[0].dataset == "/entry/data/data");
CHECK(source[0].source_first_image == 1);
CHECK(source[0].virtual_first_image == 0);
CHECK(source[0].image_count == 1);
CHECK(source[1].filename == "source_legacy_mapping_data_000002.h5");
CHECK(source[1].dataset == "/entry/data/data");
CHECK(source[1].source_first_image == 0);
CHECK(source[1].virtual_first_image == 1);
CHECK(source[1].image_count == 2);
}
remove("source_legacy_mapping_master.h5");
remove("source_legacy_mapping_data_000001.h5");
remove("source_legacy_mapping_data_000002.h5");
remove("source_legacy_mapping_data_000003.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_ProcessingHDF5_FromVDS_MapsToDataFiles", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("proc_source_vds").ImagesPerTrigger(5).ImagesPerFile(2).OverwriteExistingFiles(true);
x.BitDepthImage(16).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 51);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter writer(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
image[5678] = static_cast<int16_t>(500 + i);
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(writer.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
writer.WriteHDF5(end_message);
writer.Finalize();
}
std::vector<HDF5DataSourceMessage> source_data;
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("proc_source_vds_master.h5"));
source_data = reader.GetHDF5DataSource(1, 3);
REQUIRE(source_data.size() == 2);
CHECK(source_data[0].filename == "proc_source_vds_data_000001.h5");
CHECK(source_data[1].filename == "proc_source_vds_data_000002.h5");
}
{
DiffractionExperiment proc_x = x;
proc_x.FilePrefix("proc_from_vds")
.ImagesPerTrigger(3)
.SetFileWriterFormat(FileWriterFormat::NXmxIntegrated)
.OverwriteExistingFiles(true);
StartMessage start_message;
proc_x.FillMessage(start_message);
start_message.number_of_images = 3;
start_message.images_per_file = 3;
start_message.write_images = false;
start_message.write_master_file = true;
start_message.hdf5_source_data = source_data;
FileWriter writer(start_message);
for (int i = 0; i < 3; i++) {
DataMessage message{};
message.number = i;
message.original_number = i + 1;
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spot_count = 200 + i;
REQUIRE_NOTHROW(writer.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = 3;
writer.WriteHDF5(end_message);
writer.Finalize();
}
{
HDF5ReadOnlyFile file("proc_from_vds_master.h5");
HDF5DataSet data(file, "/entry/data/data");
HDF5Dcpl dcpl(data);
REQUIRE(dcpl.GetLayout() == HDF5DataSetLayout::VIRTUAL);
auto mappings = dcpl.GetVirtualMappings();
REQUIRE(mappings.size() == 2);
CHECK(mappings[0].filename == "proc_source_vds_data_000001.h5");
CHECK(mappings[1].filename == "proc_source_vds_data_000002.h5");
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("proc_from_vds_master.h5"));
auto img0 = reader.LoadImage(0);
REQUIRE(img0);
CHECK(img0->Image()[5678] == 501);
auto img2 = reader.LoadImage(2);
REQUIRE(img2);
CHECK(img2->Image()[5678] == 503);
}
remove("proc_source_vds_master.h5");
remove("proc_source_vds_data_000001.h5");
remove("proc_source_vds_data_000002.h5");
remove("proc_source_vds_data_000003.h5");
remove("proc_from_vds_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// rugnux --stride N processes every Nth image, so image i of the _process.h5 IS source image
// start + i*N. The mapping used to be built without the stride, linking the first N images instead,
// which put each frame's picture next to a different frame's analysis. Read the pixels back through
// the written VDS rather than only inspecting the mapping: that is what a user opening the file sees.
TEST_CASE("JFJochReader_ProcessingHDF5_Strided_LinksTheProcessedImages", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("proc_stride_src").ImagesPerTrigger(10).ImagesPerFile(2).OverwriteExistingFiles(true);
x.BitDepthImage(16).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 51);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter writer(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
image[5678] = static_cast<int16_t>(500 + i); // per-image tag: which source frame is this?
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(writer.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
writer.WriteHDF5(end_message);
writer.Finalize();
}
// start 1, 3 images, stride 3 -> source images 1, 4, 7. Two images per data file, so those sit in
// data_000001 (holds 0,1), data_000003 (holds 4,5) and data_000004 (holds 6,7).
std::vector<HDF5DataSourceMessage> source_data;
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("proc_stride_src_master.h5"));
source_data = reader.GetHDF5DataSource(1, 3, 3);
// Non-adjacent source images cannot be merged into one run, so one mapping per image.
REQUIRE(source_data.size() == 3);
CHECK(source_data[0].filename == "proc_stride_src_data_000001.h5");
CHECK(source_data[0].source_first_image == 1);
CHECK(source_data[0].virtual_first_image == 0);
CHECK(source_data[0].image_count == 1);
CHECK(source_data[1].filename == "proc_stride_src_data_000003.h5");
CHECK(source_data[1].source_first_image == 0);
CHECK(source_data[1].virtual_first_image == 1);
CHECK(source_data[1].image_count == 1);
CHECK(source_data[2].filename == "proc_stride_src_data_000004.h5");
CHECK(source_data[2].source_first_image == 1);
CHECK(source_data[2].virtual_first_image == 2);
CHECK(source_data[2].image_count == 1);
}
{
DiffractionExperiment proc_x = x;
proc_x.FilePrefix("proc_from_stride")
.ImagesPerTrigger(3)
.SetFileWriterFormat(FileWriterFormat::NXmxIntegrated)
.OverwriteExistingFiles(true);
StartMessage start_message;
proc_x.FillMessage(start_message);
start_message.number_of_images = 3;
start_message.images_per_file = 3;
start_message.write_images = false;
start_message.write_master_file = true;
start_message.hdf5_source_data = source_data;
FileWriter writer(start_message);
for (int i = 0; i < 3; i++) {
DataMessage message{};
message.number = i;
message.original_number = 1 + i * 3;
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.spot_count = 200 + i;
REQUIRE_NOTHROW(writer.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = 3;
writer.WriteHDF5(end_message);
writer.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("proc_from_stride_master.h5"));
// The whole point: frame i of the process file carries source frame 1 + 3i, not 1 + i.
auto img0 = reader.LoadImage(0);
REQUIRE(img0);
CHECK(img0->Image()[5678] == 501);
auto img1 = reader.LoadImage(1);
REQUIRE(img1);
CHECK(img1->Image()[5678] == 504);
auto img2 = reader.LoadImage(2);
REQUIRE(img2);
CHECK(img2->Image()[5678] == 507);
}
remove("proc_stride_src_master.h5");
for (int i = 1; i <= 5; i++)
remove(("proc_stride_src_data_00000" + std::to_string(i) + ".h5").c_str());
remove("proc_from_stride_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// One reflection with a distinct value in EVERY field that is meant to survive a write/read cycle,
// keyed on the image and the reflection index so no two are alike and a field read back from the
// wrong place shows up. Fractional image_number on purpose: a 3D-integrated reflection has one, and
// it is the field the offline --scale path silently lost when it was not read back.
static Reflection MakeTestReflection(int i, int j) {
const auto f = static_cast<float>(i * 10 + j);
return Reflection{
.h = 10 + i + 30 * j,
.k = 20 + j,
.l = 30 + j,
.image_number = static_cast<float>(i) + 0.25f * static_cast<float>(j + 1),
.delta_phi_deg = 0.1f + 0.01f * f,
.predicted_x = 100.0f + f,
.predicted_y = 200.0f + f,
.observed_x = 100.5f + f,
.observed_y = 200.5f + f,
.d = 1.5f + 0.1f * f,
.I = 1000.0f + f,
.bkg = 10.0f + f,
.sigma = 2.0f + 0.5f * f,
.rlp = 1.0f + 0.125f * f,
.partiality = 0.5f + 0.01f * f,
.zeta = 0.01f + 0.001f * f,
.image_scale_corr = 1.0f + 0.25f * f
};
}
// Every field of the round trip, against the reflection that was written. Deliberately NOT checked,
// because they are not part of it: dist_ewald and observed are prediction/integration scratch that
// is never written, and on_ice_ring is recomputed from the resolution by whoever scales (see the
// ice-ring handling in Rugnux and in the --scale path).
static void CheckReflectionRoundTrip(const Reflection &got, int i, int j) {
INFO("image " << i << " reflection " << j);
const Reflection want = MakeTestReflection(i, j);
CHECK(got.h == want.h);
CHECK(got.k == want.k);
CHECK(got.l == want.l);
CHECK(got.image_number == Catch::Approx(want.image_number));
CHECK(got.delta_phi_deg == Catch::Approx(want.delta_phi_deg));
CHECK(got.predicted_x == Catch::Approx(want.predicted_x));
CHECK(got.predicted_y == Catch::Approx(want.predicted_y));
CHECK(got.observed_x == Catch::Approx(want.observed_x));
CHECK(got.observed_y == Catch::Approx(want.observed_y));
CHECK(got.d == Catch::Approx(want.d));
CHECK(got.I == Catch::Approx(want.I));
CHECK(got.bkg == Catch::Approx(want.bkg));
CHECK(got.sigma == Catch::Approx(want.sigma));
CHECK(got.rlp == Catch::Approx(want.rlp)); // stored as 1/rlp, inverted again on read
CHECK(got.partiality == Catch::Approx(want.partiality));
CHECK(got.zeta == Catch::Approx(want.zeta));
CHECK(got.image_scale_corr == Catch::Approx(want.image_scale_corr));
}
TEST_CASE("JFJochReader_ReadReflections_VDS", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("read_reflections_vds")
.ImagesPerTrigger(4)
.ImagesPerFile(1)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxVDS)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 0);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter writer(start_message);
ScanResultGenerator scan_result(x);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
if (i == 1 || i == 3) {
message.integrated_reflections = 2;
message.reflections = {MakeTestReflection(i, 0), MakeTestReflection(i, 1)};
message.mosaicity_deg = i*0.15f;
message.indexing_lattice = CrystalLattice({100,0,0}, {0,50,0}, {0,0,30});
}
REQUIRE_NOTHROW(writer.WriteHDF5(message));
scan_result.Add(message);
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
scan_result.FillEndMessage(end_message);
writer.WriteHDF5(end_message);
writer.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_reflections_vds_master.h5"));
auto reflections = reader.ReadReflections();
REQUIRE(reflections.size() == 4);
CHECK(reflections[0].reflections.empty());
REQUIRE(reflections[1].reflections.size() == 2);
CheckReflectionRoundTrip(reflections[1].reflections[0], 1, 0);
CheckReflectionRoundTrip(reflections[1].reflections[1], 1, 1);
CHECK(reflections[1].mosaicity_deg == Catch::Approx(0.15f));
CHECK(reflections[1].latt.CalcVolume() == Catch::Approx(100*50*30));
CHECK(reflections[2].reflections.empty());
REQUIRE(reflections[3].reflections.size() == 2);
CheckReflectionRoundTrip(reflections[3].reflections[0], 3, 0);
CheckReflectionRoundTrip(reflections[3].reflections[1], 3, 1);
CHECK(reflections[3].mosaicity_deg == Catch::Approx(0.45f));
CHECK(reflections[3].latt.Vec0().x == Catch::Approx(100.0f));
CHECK(reflections[3].latt.Vec1().y == Catch::Approx(50.0f));
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_reflections_vds_master.h5"));
auto reflections = reader.ReadReflections(1, 3);
REQUIRE(reflections.size() == 3);
REQUIRE(reflections[0].reflections.size() == 2); // original image 1
CheckReflectionRoundTrip(reflections[0].reflections[0], 1, 0);
CHECK(reflections[1].reflections.empty()); // original image 2
REQUIRE(reflections[2].reflections.size() == 2); // original image 3
CheckReflectionRoundTrip(reflections[2].reflections[0], 3, 0);
}
remove("read_reflections_vds_master.h5");
remove("read_reflections_vds_data_000001.h5");
remove("read_reflections_vds_data_000002.h5");
remove("read_reflections_vds_data_000003.h5");
remove("read_reflections_vds_data_000004.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// The per-image reflections and lattices are written in the setting the images were indexed in, but
// the space group - and with it the conventional setting the cell beside them is in - is only settled
// after the merge, so the two can differ by an integral change of basis. /entry/MX/reindexMatrix
// carries it, and the reader applies it, so what comes out is in the cell's setting. A file without
// the dataset (every file written before it existed) is read as the identity, which is what the
// round-trip tests above check.
TEST_CASE("JFJochReader_ReadReflections_Reindex", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("read_reflections_reindex")
.ImagesPerTrigger(2)
.ImagesPerFile(1)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxVDS)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
// hkl_cell = M . hkl_written, det 2 - the size of step a primitive-to-centred re-seat takes.
const std::array<int32_t, 9> M = {1, 1, 0,
0, 1, 1,
1, 0, 1};
std::vector<int16_t> image(x.GetPixelsNum(), 0);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter writer(start_message);
ScanResultGenerator scan_result(x);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
if (i == 1) {
message.integrated_reflections = 2;
message.reflections = {MakeTestReflection(i, 0), MakeTestReflection(i, 1)};
message.indexing_result = true;
message.indexing_lattice = CrystalLattice({100,0,0}, {0,50,0}, {0,0,30});
}
REQUIRE_NOTHROW(writer.WriteHDF5(message));
scan_result.Add(message);
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
end_message.reindex_matrix = M;
scan_result.FillEndMessage(end_message);
writer.WriteHDF5(end_message);
writer.Finalize();
}
// hkl and the lattice come back in the cell's setting; every other field is untouched.
const auto check = [&](const Reflection &got, int j) {
const Reflection want = MakeTestReflection(1, j);
INFO("reflection " << j);
CHECK(got.h == M[0] * want.h + M[1] * want.k + M[2] * want.l);
CHECK(got.k == M[3] * want.h + M[4] * want.k + M[5] * want.l);
CHECK(got.l == M[6] * want.h + M[7] * want.k + M[8] * want.l);
CHECK(got.I == Catch::Approx(want.I));
CHECK(got.d == Catch::Approx(want.d));
CHECK(got.image_number == Catch::Approx(want.image_number));
};
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_reflections_reindex_master.h5"));
REQUIRE(reader.GetDataset()->reindex_matrix.has_value());
CHECK(reader.GetDataset()->reindex_matrix.value() == M);
auto reflections = reader.ReadReflections();
REQUIRE(reflections.size() == 2);
REQUIRE(reflections[1].reflections.size() == 2);
check(reflections[1].reflections[0], 0);
check(reflections[1].reflections[1], 1);
// latt = M . latt_written, row by row: (100,0,0)+(0,50,0), (0,50,0)+(0,0,30), (100,0,0)+(0,0,30).
CHECK(reflections[1].latt.Vec0().x == Catch::Approx(100.0f));
CHECK(reflections[1].latt.Vec0().y == Catch::Approx(50.0f));
CHECK(reflections[1].latt.Vec1().y == Catch::Approx(50.0f));
CHECK(reflections[1].latt.Vec1().z == Catch::Approx(30.0f));
CHECK(reflections[1].latt.Vec2().x == Catch::Approx(100.0f));
CHECK(reflections[1].latt.Vec2().z == Catch::Approx(30.0f));
CHECK(reflections[1].latt.CalcVolume() == Catch::Approx(2.0 * 100 * 50 * 30));
// The per-image message path (the viewer's) is re-seated the same way.
auto reader_image = reader.LoadImage(1);
REQUIRE(reader_image);
REQUIRE(reader_image->ImageData().reflections.size() == 2);
check(reader_image->ImageData().reflections[0], 0);
REQUIRE(reader_image->ImageData().indexing_lattice);
CHECK(reader_image->ImageData().indexing_lattice->CalcVolume()
== Catch::Approx(2.0 * 100 * 50 * 30));
}
remove("read_reflections_reindex_master.h5");
remove("read_reflections_reindex_data_000001.h5");
remove("read_reflections_reindex_data_000002.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
static std::vector<SpotToSave> MakeTestSpots(int i) {
return {
SpotToSave{
.x = 1, .y = 2, .intensity = 376,
.h = 11, .k = -3, .l = -5,
.dist_ewald_sphere = 0.1234f,
.ice_ring = true,
.indexed = true
},
SpotToSave{
.x = 7, .y = static_cast<float>(-3 - i), .intensity = 0.156f,
.ice_ring = false,
.indexed = false,
}
};
}
// Assert the full field set on spots[0] and the per-image variation on
// spots[1].y, which is the only field that differs across images.
static void CheckSpotFields(const SpotToSave &s0, const SpotToSave &s1, int i) {
CHECK(s0.x == 1);
CHECK(s0.y == 2);
CHECK(s0.intensity == Catch::Approx(376));
CHECK(s0.ice_ring == true);
CHECK(s0.indexed == true);
CHECK(s0.h == 11);
CHECK(s0.k == -3);
CHECK(s0.l == -5);
CHECK(s0.dist_ewald_sphere == Catch::Approx(0.1234f));
CHECK(s0.image == i);
CHECK(s1.x == Catch::Approx(7));
CHECK(s1.y == Catch::Approx(static_cast<float>(-3 - i)));
CHECK(s1.intensity == Catch::Approx(0.156f));
CHECK(s1.ice_ring == false);
CHECK(s1.indexed == false);
CHECK(s1.image == i);
}
TEST_CASE("JFJochReader_ReadSpots_Legacy", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("read_spots_legacy")
.ImagesPerTrigger(4)
.ImagesPerFile(1)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxLegacy)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 0);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
message.spots = MakeTestSpots(i);
message.spot_count = 72 + i;
message.spot_count_ice_rings = 45 + 2 * i;
message.spot_count_low_res = 12 + 3 * i;
message.spot_count_indexed = 15 + 4 * i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
// All images, one at a time
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_legacy_master.h5"));
for (int i = 0; i < 4; i++) {
std::vector<SpotToSave> spots;
REQUIRE_NOTHROW(spots = reader.ReadSpots(i));
REQUIRE(spots.size() == 2);
CheckSpotFields(spots[0], spots[1], i);
}
}
// Out-of-range must throw
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_legacy_master.h5"));
REQUIRE_THROWS(reader.ReadSpots(4));
}
remove("read_spots_legacy_master.h5");
remove("read_spots_legacy_data_000001.h5");
remove("read_spots_legacy_data_000002.h5");
remove("read_spots_legacy_data_000003.h5");
remove("read_spots_legacy_data_000004.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_ReadSpots_VDS", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("read_spots_vds")
.ImagesPerTrigger(4)
.ImagesPerFile(1)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxVDS)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 0);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
message.spots = MakeTestSpots(i);
message.spot_count = 72 + i;
message.spot_count_ice_rings = 45 + 2 * i;
message.spot_count_low_res = 12 + 3 * i;
message.spot_count_indexed = 15 + 4 * i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
// All images, one at a time — also verifies that .image carries the
// correct global index across the virtual-to-source remapping.
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_vds_master.h5"));
for (int i = 0; i < 4; i++) {
std::vector<SpotToSave> spots;
REQUIRE_NOTHROW(spots = reader.ReadSpots(i));
REQUIRE(spots.size() == 2);
CheckSpotFields(spots[0], spots[1], i);
}
}
// Out-of-range must throw
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_vds_master.h5"));
REQUIRE_THROWS(reader.ReadSpots(4));
}
// Image with no spots returns an empty vector
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_vds_master.h5"));
// Write a separate 2-image VDS file where only image 0 has spots.
DiffractionExperiment y(DetJF(1));
y.FilePrefix("read_spots_vds_sparse")
.ImagesPerTrigger(2)
.ImagesPerFile(1)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxVDS)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
{
StartMessage start_message;
y.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < 2; i++) {
DataMessage message{};
message.image = CompressedImage(image, y.GetXPixelsNum(), y.GetYPixelsNum());
message.number = i;
if (i == 0)
message.spots = MakeTestSpots(i);
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = 2;
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
JFJochHDF5Reader sparse_reader;
REQUIRE_NOTHROW(sparse_reader.ReadFile("read_spots_vds_sparse_master.h5"));
std::vector<SpotToSave> spots_0, spots_1;
REQUIRE_NOTHROW(spots_0 = sparse_reader.ReadSpots(0));
REQUIRE_NOTHROW(spots_1 = sparse_reader.ReadSpots(1));
REQUIRE(spots_0.size() == 2);
CHECK(spots_0[0].image == 0);
CHECK(spots_1.empty());
remove("read_spots_vds_sparse_master.h5");
remove("read_spots_vds_sparse_data_000001.h5");
remove("read_spots_vds_sparse_data_000002.h5");
}
remove("read_spots_vds_master.h5");
remove("read_spots_vds_data_000001.h5");
remove("read_spots_vds_data_000002.h5");
remove("read_spots_vds_data_000003.h5");
remove("read_spots_vds_data_000004.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_ReadAllSpots_VDS", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("read_spots_vds")
.ImagesPerTrigger(20)
.ImagesPerFile(3)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxVDS)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 0);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
message.spots = MakeTestSpots(i);
message.spot_count = 72 + i;
message.spot_count_ice_rings = 45 + 2 * i;
message.spot_count_low_res = 12 + 3 * i;
message.spot_count_indexed = 15 + 4 * i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_vds_master.h5"));
std::shared_ptr<JFJochReaderSpots> ret;
REQUIRE_NOTHROW(ret = reader.ReadAllSpots(1, 15, 2));
// 1,3,5,7,9,11,13,15
REQUIRE(ret);
REQUIRE(ret->start_image == 1);
REQUIRE(ret->stride == 2);
REQUIRE(ret->spots.size() == 8);
for (int i = 0; i < ret->spots.size(); i++) {
REQUIRE(ret->spots[i].size() == 2);
CheckSpotFields(ret->spots[i][0], ret->spots[i][1], 2 * i + 1);
}
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_vds_master.h5"));
REQUIRE_THROWS(reader.ReadAllSpots(-5,0));
REQUIRE_THROWS(reader.ReadAllSpots(5,0));
}
remove("read_spots_vds_master.h5");
remove("read_spots_vds_data_000001.h5");
remove("read_spots_vds_data_000002.h5");
remove("read_spots_vds_data_000003.h5");
remove("read_spots_vds_data_000004.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_ReadSpots_Integrated", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.FilePrefix("read_spots_integrated")
.ImagesPerTrigger(4)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxIntegrated)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
std::vector<int16_t> image(x.GetPixelsNum(), 0);
RegisterHDF5Filter();
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < x.GetImageNum(); i++) {
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
message.spots = MakeTestSpots(i);
message.spot_count = 72 + i;
message.spot_count_ice_rings = 45 + 2 * i;
message.spot_count_low_res = 12 + 3 * i;
message.spot_count_indexed = 15 + 4 * i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
// All images, one at a time
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_integrated_master.h5"));
for (int i = 0; i < 4; i++) {
std::vector<SpotToSave> spots;
REQUIRE_NOTHROW(spots = reader.ReadSpots(i));
REQUIRE(spots.size() == 2);
CheckSpotFields(spots[0], spots[1], i);
}
}
// Out-of-range must throw
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_integrated_master.h5"));
REQUIRE_THROWS(reader.ReadSpots(4));
}
// Image with no spots returns an empty vector
{
DiffractionExperiment y(DetJF(1));
y.FilePrefix("read_spots_integrated_sparse")
.ImagesPerTrigger(3)
.OverwriteExistingFiles(true)
.BitDepthImage(16)
.PixelSigned(true)
.SetFileWriterFormat(FileWriterFormat::NXmxIntegrated)
.IndexingAlgorithm(IndexingAlgorithmEnum::FFT)
.Compression(CompressionAlgorithm::NO_COMPRESSION);
{
StartMessage start_message;
y.FillMessage(start_message);
FileWriter file_set(start_message);
for (int i = 0; i < 3; i++) {
DataMessage message{};
message.image = CompressedImage(image, y.GetXPixelsNum(), y.GetYPixelsNum());
message.number = i;
if (i == 1)
message.spots = MakeTestSpots(i);
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = 3;
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("read_spots_integrated_sparse_master.h5"));
CHECK(reader.ReadSpots(0).empty());
REQUIRE(reader.ReadSpots(1).size() == 2);
CheckSpotFields(reader.ReadSpots(1)[0], reader.ReadSpots(1)[1], 1);
CHECK(reader.ReadSpots(2).empty());
remove("read_spots_integrated_sparse_master.h5");
}
remove("read_spots_integrated_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
TEST_CASE("JFJochReader_Snapshots", "[HDF5][Full]") {
RegisterHDF5Filter();
DiffractionExperiment x(DetJF(1));
x.FilePrefix("test_snap").ImagesPerTrigger(4).OverwriteExistingFiles(true);
x.BitDepthImage(16).ImagesPerFile(1).SetFileWriterFormat(FileWriterFormat::NXmxVDS).PixelSigned(true);
x.Compression(CompressionAlgorithm::NO_COMPRESSION);
// 1. Original dataset: distinct pixels + "original" MX metadata (not indexed, bkg = 10 + i).
std::vector<int16_t> image(x.GetPixelsNum());
{
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
ScanResultGenerator generator(x);
for (int i = 0; i < 4; i++) {
image[5678] = 100 + i;
DataMessage message{};
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
message.indexing_result = false;
message.bkg_estimate = 10.0 + i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
}
EndMessage end_message;
end_message.max_image_number = 4;
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
// 2. A reprocessing result over the same 4 images: integrated master, "reprocessed" MX
// metadata (all indexed, bkg = 99 + i) and deliberately wrong pixels that must NOT surface
// (snapshot pixels still come from the original image source).
{
DiffractionExperiment px(x);
px.FilePrefix("test_snap_proc").SetFileWriterFormat(FileWriterFormat::NXmxIntegrated);
StartMessage start_message;
px.FillMessage(start_message);
FileWriter file_set(start_message);
ScanResultGenerator generator(px);
std::vector<int16_t> proc_image(x.GetPixelsNum(), 7);
for (int i = 0; i < 4; i++) {
DataMessage message{};
message.image = CompressedImage(proc_image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
message.indexing_result = true;
message.bkg_estimate = 99.0 + i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
generator.Add(message);
}
EndMessage end_message;
end_message.max_image_number = 4;
generator.FillEndMessage(end_message);
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_snap_master.h5"));
REQUIRE(reader.GetNumberOfImages() == 4);
CHECK(reader.ActiveSnapshot() == "Original");
CHECK(reader.SnapshotNames() == std::vector<std::string>{"Original"});
// Dataset-level plot arrays come from the active (original) metadata source.
REQUIRE(reader.GetDataset()->bkg_estimate.size() == 4);
CHECK(reader.GetDataset()->bkg_estimate[2] == Catch::Approx(12.0));
auto orig0 = reader.LoadImage(0);
REQUIRE(orig0);
CHECK(orig0->Image()[5678] == 100);
REQUIRE(orig0->ImageData().indexing_result.has_value());
CHECK(orig0->ImageData().indexing_result.value() == false);
CHECK(orig0->ImageData().bkg_estimate.value() == Catch::Approx(10.0));
// Register the reprocessing result as a second metadata source over the same images.
REQUIRE_NOTHROW(reader.RegisterSnapshot("Reprocess", "test_snap_proc_master.h5"));
{
auto names = reader.SnapshotNames();
CHECK(std::find(names.begin(), names.end(), "Original") != names.end());
CHECK(std::find(names.begin(), names.end(), "Reprocess") != names.end());
}
REQUIRE_NOTHROW(reader.SetActiveSnapshot("Reprocess"));
CHECK(reader.ActiveSnapshot() == "Reprocess");
// Plots now come from the reprocessing master.
CHECK(reader.GetDataset()->bkg_estimate[2] == Catch::Approx(101.0));
auto repro0 = reader.LoadImage(0);
REQUIRE(repro0);
// Pixels still from the original image source, not the 7's stored in the process file.
CHECK(repro0->Image()[5678] == 100);
// Metadata from the reprocessing snapshot.
CHECK(repro0->ImageData().indexing_result.value() == true);
CHECK(repro0->ImageData().bkg_estimate.value() == Catch::Approx(99.0));
// Switch back to the original metadata.
REQUIRE_NOTHROW(reader.SetActiveSnapshot("Original"));
auto orig0b = reader.LoadImage(0);
REQUIRE(orig0b);
CHECK(orig0b->Image()[5678] == 100);
CHECK(orig0b->ImageData().indexing_result.value() == false);
REQUIRE_THROWS(reader.SetActiveSnapshot("Nonexistent"));
}
remove("test_snap_master.h5");
remove("test_snap_data_000001.h5");
remove("test_snap_data_000002.h5");
remove("test_snap_data_000003.h5");
remove("test_snap_data_000004.h5");
remove("test_snap_proc_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// The chain may be sent in the END message or left to the writer to build. Both must produce the
// same file - the sent one is written verbatim, which is what will later allow measured positions to
// be reported, and the built one is what a producer that does not send it gets.
TEST_CASE("JFJochReader_TransformationChain_SentAndBuilt", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(5).OverwriteExistingFiles(true);
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
x.Goniometer(GoniometerAxis("omega", 95, 0.1f, Coord(0,-1,0), {}));
x.Smargon(SmargonPosition{.phi_deg = -7.25f, .chi_deg = 12.5f});
RegisterHDF5Filter();
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
const auto write = [&](const std::string &prefix, bool send_chain) {
DiffractionExperiment local = x;
local.FilePrefix(prefix);
StartMessage start_message;
local.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < local.GetImageNum(); i++) {
message.image = CompressedImage(image, local.GetXPixelsNum(), local.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = local.GetImageNum();
if (send_chain)
end_message.transformations = local.BuildTransformationChain(local.GetImageNum());
file_set.WriteHDF5(end_message);
file_set.Finalize();
};
write("test_chain_built", false);
write("test_chain_sent", true);
const auto read = [](const std::string &prefix) {
JFJochHDF5Reader reader;
reader.ReadFile(prefix + "_master.h5");
return reader.GetDataset()->experiment;
};
const auto built = read("test_chain_built");
const auto sent = read("test_chain_sent");
REQUIRE(built.GetGoniometer().has_value());
REQUIRE(sent.GetGoniometer().has_value());
CHECK(sent.GetGoniometer()->GetName() == built.GetGoniometer()->GetName());
CHECK(sent.GetGoniometer()->GetStart_deg()
== Catch::Approx(built.GetGoniometer()->GetStart_deg()).margin(1e-3));
CHECK(sent.GetGoniometer()->GetIncrement_deg()
== Catch::Approx(built.GetGoniometer()->GetIncrement_deg()).margin(1e-4));
// chi/phi survive both routes, which they did not before they became ordinary axes.
REQUIRE(built.GetDatasetSettings().GetSmargonPosition().has_value());
REQUIRE(sent.GetDatasetSettings().GetSmargonPosition().has_value());
CHECK(sent.GetDatasetSettings().GetSmargonPosition()->chi_deg
== Catch::Approx(12.5f).margin(1e-3));
CHECK(sent.GetDatasetSettings().GetSmargonPosition()->phi_deg
== Catch::Approx(-7.25f).margin(1e-3));
// Compared on the files, not through the reader: the reader reads neither AXISNAME_end nor the
// rotation width, so it cannot see the two routes diverge - and it did, until the writer started
// deriving them for a chain it was handed.
{
HDF5ReadOnlyFile built_file("test_chain_built_master.h5");
HDF5ReadOnlyFile sent_file("test_chain_sent_master.h5");
CHECK(built_file.FindLeafs("/entry/sample/transformations")
== sent_file.FindLeafs("/entry/sample/transformations"));
CHECK(sent_file.ReadVector<double>("/entry/sample/transformations/omega_end").size()
== static_cast<size_t>(x.GetImageNum()));
CHECK(sent_file.ReadVector<double>("/entry/sample/transformations/omega_range_average").at(0)
== Catch::Approx(0.1).margin(1e-4));
}
remove("test_chain_built_master.h5");
remove("test_chain_sent_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// Recovering the sample axes from a written file, across the configurations the writer produces.
// The reader searches every leaf of /entry/sample/transformations, so on the way it meets the
// writer's own AXISNAME_end and rotation-width datasets, which are not axes and carry no
// transformation_type. It used to throw on them: a master whose axis did not turn never stopped the
// search early, walked into omega_end and could not be opened at all - which took out rugnux's own
// output for a grid scan.
TEST_CASE("JFJochReader_AxisRecovery", "[HDF5][Full]") {
RegisterHDF5Filter();
const auto round_trip = [](DiffractionExperiment x, const std::string &prefix) {
x.FilePrefix(prefix).OverwriteExistingFiles(true)
.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
StartMessage start_message;
x.FillMessage(start_message);
{
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < x.GetImageNum(); i++) {
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
}
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile(prefix + "_master.h5"));
return reader.GetDataset()->experiment;
};
const auto cleanup = [](const std::string &prefix) {
remove((prefix + "_master.h5").c_str());
remove((prefix + "_data_000001.h5").c_str());
};
SECTION("a sweep") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(5).Goniometer(GoniometerAxis("omega", 95, 0.1f, Coord(0,-1,0), {}));
const auto out = round_trip(x, "test_ax_sweep");
REQUIRE(out.GetGoniometer().has_value());
CHECK(out.GetGoniometer()->GetName() == "omega");
CHECK(out.GetGoniometer()->IsScanning());
CHECK(out.GetGoniometer()->GetStart_deg() == Catch::Approx(95).margin(1e-3));
CHECK(out.GetGoniometer()->GetIncrement_deg() == Catch::Approx(0.1).margin(1e-4));
CHECK(out.GetGoniometer()->GetAxis() == Coord(0,-1,0));
cleanup("test_ax_sweep");
}
SECTION("a sweep about an axis that is not called omega") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(5).Goniometer(GoniometerAxis("kappa", 10, 0.5f, Coord(0,-1,0), {}));
const auto out = round_trip(x, "test_ax_kappa");
REQUIRE(out.GetGoniometer().has_value());
CHECK(out.GetGoniometer()->GetName() == "kappa");
CHECK(out.GetGoniometer()->IsScanning());
CHECK(out.GetGoniometer()->GetIncrement_deg() == Catch::Approx(0.5).margin(1e-4));
cleanup("test_ax_kappa");
}
SECTION("a spindle that does not turn") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(5).Goniometer(GoniometerAxis("omega", 12.5f, 0.0f, Coord(0,-1,0), {}));
const auto out = round_trip(x, "test_ax_still");
REQUIRE(out.GetGoniometer().has_value());
CHECK(out.GetGoniometer()->GetName() == "omega");
CHECK(!out.GetGoniometer()->IsScanning());
CHECK(out.GetGoniometer()->GetStart_deg() == Catch::Approx(12.5).margin(1e-3));
cleanup("test_ax_still");
}
SECTION("a grid scan, which sits on a spindle that does not turn") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(6).GridScan(GridScanSettings(3, 10.0f, 20.0f, false, false).ImageNum(6));
const auto out = round_trip(x, "test_ax_grid");
REQUIRE(out.GetGridScan().has_value());
CHECK(out.GetGridScan()->GetNFast() == 3);
CHECK(out.GetGridScan()->GetGridStepX_um() == Catch::Approx(10.0).margin(1e-3));
CHECK(out.GetGridScan()->GetGridStepY_um() == Catch::Approx(20.0).margin(1e-3));
// NXmx cannot say "no rotation", so the writer records the spindle standing still.
REQUIRE(out.GetGoniometer().has_value());
CHECK(!out.GetGoniometer()->IsScanning());
cleanup("test_ax_grid");
}
SECTION("a grid scan under a turning spindle") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(6).GridScan(GridScanSettings(3, 10.0f, 20.0f, false, false).ImageNum(6))
.Goniometer(GoniometerAxis("omega", 0, 0.2f, Coord(0,-1,0), {}));
const auto out = round_trip(x, "test_ax_gridsweep");
REQUIRE(out.GetGridScan().has_value());
CHECK(out.GetGridScan()->GetNFast() == 3);
REQUIRE(out.GetGoniometer().has_value());
CHECK(out.GetGoniometer()->IsScanning());
CHECK(out.GetGoniometer()->GetIncrement_deg() == Catch::Approx(0.2).margin(1e-4));
cleanup("test_ax_gridsweep");
}
SECTION("a sweep with the head at a Smargon position") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(5).Goniometer(GoniometerAxis("omega", 95, 0.1f, Coord(0,-1,0), {}))
.Smargon(SmargonPosition{.phi_deg = -7.25f, .chi_deg = 12.5f});
const auto out = round_trip(x, "test_ax_smargon");
REQUIRE(out.GetGoniometer().has_value());
CHECK(out.GetGoniometer()->GetName() == "omega");
CHECK(out.GetGoniometer()->IsScanning());
REQUIRE(out.GetDatasetSettings().GetSmargonPosition().has_value());
CHECK(out.GetDatasetSettings().GetSmargonPosition()->chi_deg == Catch::Approx(12.5).margin(1e-3));
CHECK(out.GetDatasetSettings().GetSmargonPosition()->phi_deg == Catch::Approx(-7.25).margin(1e-3));
cleanup("test_ax_smargon");
}
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// saturation_value is written inclusive and used exclusive, so a read has to add the count back. It
// did not, and the value fell by one on every write-read-write cycle - unbounded, and compounding
// whenever a _process.h5 was reprocessed. Nothing caught it: no test asserted the read-back limit.
TEST_CASE("JFJochReader_SaturationSurvivesRoundTrip", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(2).OverwriteExistingFiles(true).FilePrefix("test_satrt");
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
const int64_t original_limit = x.GetSaturationLimit();
RegisterHDF5Filter();
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
// The second pass writes metadata only: reading pins the experiment to signed 32-bit, the
// container images are handed out in, so feeding it the uint16 frames again would - rightly -
// be refused by the writer's pixel-format check.
const auto write = [&](const DiffractionExperiment &src, const std::string &prefix,
bool with_images) {
DiffractionExperiment local = src;
local.FilePrefix(prefix).OverwriteExistingFiles(true);
StartMessage start_message;
local.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
if (with_images) {
for (int i = 0; i < 2; i++) {
message.image = CompressedImage(image, local.GetXPixelsNum(), local.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
}
EndMessage end_message;
end_message.max_image_number = 2;
file_set.WriteHDF5(end_message);
file_set.Finalize();
return start_message.saturation_value;
};
const auto read = [](const std::string &prefix) {
JFJochHDF5Reader reader;
reader.ReadFile(prefix + "_master.h5");
return reader.GetDataset()->experiment;
};
const int64_t declared_once = write(x, "test_satrt", true);
CHECK(declared_once == SaturationValueFromLimit(original_limit));
const auto once = read("test_satrt");
CHECK(once.GetSaturationLimit() == original_limit);
// The cycle that used to lose a count: read a file, write what was read, read it again.
const int64_t declared_twice = write(once, "test_satrt2", false);
CHECK(declared_twice == declared_once);
CHECK(read("test_satrt2").GetSaturationLimit() == original_limit);
remove("test_satrt_master.h5");
remove("test_satrt_data_000001.h5");
remove("test_satrt2_master.h5");
remove("test_satrt2_data_000001.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// A head position is not a sweep, and it is not the spindle either. Both properties are carried by
// the file itself - the axis length says how many images there are, the equipment_component tag says
// what the axis is - so both are checked here on the file, not through the reader: the reader alone
// cannot see a shape or an attribute it never looks at.
TEST_CASE("JFJochReader_Smargon_StillIsNotOneImage", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(5).OverwriteExistingFiles(true).FilePrefix("test_smargon");
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
x.Smargon(SmargonPosition{.phi_deg = -7.25f, .chi_deg = 12.5f});
REQUIRE(!x.GetGoniometer().has_value());
RegisterHDF5Filter();
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < x.GetImageNum(); i++) {
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
{
HDF5ReadOnlyFile master("test_smargon_master.h5");
// One entry per image. A reader takes the image count from the innermost axis of the sample
// chain when no axis varies; as scalars these read back as a single image.
CHECK(master.GetDimension("/entry/sample/transformations/chi")
== std::vector<hsize_t>{static_cast<hsize_t>(x.GetImageNum())});
CHECK(master.GetDimension("/entry/sample/transformations/phi")
== std::vector<hsize_t>{static_cast<hsize_t>(x.GetImageNum())});
CHECK(master.ReadVector<double>("/entry/sample/transformations/phi")
== std::vector<double>(x.GetImageNum(), -7.25));
// Tagged, so neither is mistaken for the spindle - and so a phi from anywhere else is not
// mistaken for a head position.
HDF5DataSet chi(master, "/entry/sample/transformations/chi");
HDF5DataSet phi(master, "/entry/sample/transformations/phi");
REQUIRE(chi.AttrExists("equipment_component"));
REQUIRE(phi.AttrExists("equipment_component"));
CHECK(chi.ReadAttrStr("equipment_component") == "smargon");
CHECK(phi.ReadAttrStr("equipment_component") == "smargon");
}
const auto read = [](const std::string &prefix) {
JFJochHDF5Reader reader;
reader.ReadFile(prefix + "_master.h5");
return reader.GetDataset()->experiment;
};
const auto read_back = read("test_smargon");
// chi is the alphabetically first stationary axis in the file; it must not become the spindle.
CHECK(!read_back.GetGoniometer().has_value());
REQUIRE(read_back.GetDatasetSettings().GetSmargonPosition().has_value());
CHECK(read_back.GetDatasetSettings().GetSmargonPosition()->chi_deg == Catch::Approx(12.5f).margin(1e-3));
CHECK(read_back.GetDatasetSettings().GetSmargonPosition()->phi_deg == Catch::Approx(-7.25f).margin(1e-3));
remove("test_smargon_master.h5");
remove("test_smargon_data_000001.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// phi is an ordinary spindle name in MX. A file whose rotation axis is called phi carries no
// equipment_component, so it stays the spindle and no head position is invented from it - which also
// means the file can be written back out, instead of colliding on a second dataset called phi.
TEST_CASE("JFJochReader_Goniometer_NamedPhiIsNotSmargon", "[HDF5][Full]") {
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(5).OverwriteExistingFiles(true).FilePrefix("test_phispindle");
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
x.Goniometer(GoniometerAxis("phi", 30, 0.2f, Coord(0,-1,0), {}));
RegisterHDF5Filter();
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < x.GetImageNum(); i++) {
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
const auto read = [](const std::string &prefix) {
JFJochHDF5Reader reader;
reader.ReadFile(prefix + "_master.h5");
return reader.GetDataset()->experiment;
};
const auto read_back = read("test_phispindle");
REQUIRE(read_back.GetGoniometer().has_value());
CHECK(read_back.GetGoniometer()->GetName() == "phi");
CHECK(read_back.GetGoniometer()->GetStart_deg() == Catch::Approx(30).margin(1e-3));
CHECK(!read_back.GetDatasetSettings().GetSmargonPosition().has_value());
// Writing what was read must not try to create phi a second time.
DiffractionExperiment rewrite = read_back;
rewrite.FilePrefix("test_phispindle_out").OverwriteExistingFiles(true);
StartMessage out_start;
rewrite.FillMessage(out_start);
FileWriter out(out_start);
EndMessage out_end;
out_end.max_image_number = rewrite.GetImageNum();
REQUIRE_NOTHROW(out.WriteHDF5(out_end));
REQUIRE_NOTHROW(out.Finalize());
remove("test_phispindle_master.h5");
remove("test_phispindle_data_000001.h5");
remove("test_phispindle_out_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// A valid NXmx master written outside the DECTRIS toolchain, in the shape a Diamond-written one
// takes: lengths in millimetres, no detectorSpecific, the distance one level up in NXinstrument, a
// pixel_mask that is an external link into a file not holding it, and per-file links naming a plain
// /data rather than /entry/data/data. Every one of those was enough on its own to stop the file
// opening, and the last one did it without an error - the run reported no images and succeeded.
namespace {
void WriteThirdPartyDataFile(const std::string &filename, const std::vector<uint16_t> &image,
hsize_t nimages, hsize_t ny, hsize_t nx) {
std::vector<uint16_t> block;
for (hsize_t i = 0; i < nimages; i++)
block.insert(block.end(), image.begin(), image.end());
HDF5File file(filename);
file.SaveVector("/data", block, {nimages, ny, nx});
}
}
TEST_CASE("JFJochReader_ThirdPartyNXmxMaster", "[HDF5][Full]") {
const hsize_t nx = 8, ny = 6, per_file = 2;
std::vector<uint16_t> image(nx * ny);
for (size_t i = 0; i < image.size(); i++)
image[i] = static_cast<uint16_t>(i * 3 + 1);
WriteThirdPartyDataFile("third_party_000001.h5", image, per_file, ny, nx);
WriteThirdPartyDataFile("third_party_000002.h5", image, per_file, ny, nx);
{
HDF5File master("third_party_master.h5");
HDF5Group entry(master, "entry");
entry.SaveScalar("definition", "NXmx");
HDF5Group instrument(entry, "instrument");
// The distance NXdetector does not carry, in millimetres
instrument.SaveScalar("detector_distance", 287.5)->Units("mm");
HDF5Group beam(instrument, "beam");
beam.SaveScalar("incident_wavelength", 0.9794)->Units("angstrom");
HDF5Group detector(instrument, "detector");
detector.SaveScalar("description", "Eiger 16M");
detector.SaveScalar("beam_center_x", 4.0)->Units("pixels");
detector.SaveScalar("beam_center_y", 3.0)->Units("pixels");
detector.SaveScalar("count_time", 0.01);
detector.SaveScalar("saturation_value", static_cast<int64_t>(65535));
detector.SaveScalar("x_pixel_size", 0.075)->Units("mm");
detector.SaveScalar("y_pixel_size", 0.075)->Units("mm");
detector.SaveScalar("sensor_thickness", 0.45)->Units("mm");
// Links into a file that does not exist at all, so neither can be dereferenced
detector.ExternalLink("third_party_no_such_meta.h5", "/mask", "pixel_mask");
HDF5Group data(entry, "data");
data.ExternalLink("third_party_000001.h5", "/data", "data_000001");
data.ExternalLink("third_party_000002.h5", "/data", "data_000002");
}
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("third_party_master.h5"));
auto dataset = reader.GetDataset();
// Images found through the link's own target path, and sized from the array itself
CHECK(dataset->experiment.GetImageNum() == 2 * per_file);
CHECK(dataset->experiment.GetXPixelsNum() == nx);
CHECK(dataset->experiment.GetYPixelsNum() == ny);
// Millimetres read as millimetres
CHECK(dataset->experiment.GetDetectorDistance_mm() == Catch::Approx(287.5));
CHECK(dataset->experiment.GetDetectorSetup().GetPixelSize_mm() == Catch::Approx(0.075));
CHECK(dataset->experiment.GetDetectorSetup().GetSensorThickness_um() == Catch::Approx(450.0));
// Both mask links dangle; the reader must fall back to an empty mask, not throw
REQUIRE(dataset->pixel_mask);
std::shared_ptr<JFJochReaderRawImage> reader_image;
for (int i = 0; i < 2 * static_cast<int>(per_file); i++) {
REQUIRE_NOTHROW(reader_image = reader.GetRawImage(i));
CHECK(reader_image->image.GetWidth() == nx);
CHECK(reader_image->image.GetHeight() == ny);
}
}
// A master that names data files nothing can be read from is an error, not an empty data set
remove("third_party_000001.h5");
remove("third_party_000002.h5");
{
JFJochHDF5Reader reader;
REQUIRE_THROWS(reader.ReadFile("third_party_master.h5"));
}
remove("third_party_master.h5");
// No leftover HDF5 objects
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// A detector swung out on a 2theta arm. NXmx has no field for it: the swing is one rotation in the
// depends_on chain the detector's position is stated as, and "two_theta" is only one beamline's name
// for that dataset. So the chain is what the reader follows, and the chain here carries two rotations
// about different axes, outboard of the translation that sets the distance - a file that stated only
// the innermost one, or composed them the other way round, gives a different plane.
//
// Both axes are stated in McStas, which is the internal frame turned half a turn about z: a reader
// that takes the vector as it stands swings the detector the wrong way, which is twice the error of
// not reading it at all.
TEST_CASE("JFJochReader_DetectorTwoThetaArm", "[HDF5][Full]") {
const hsize_t nx = 8, ny = 6;
const double two_theta_deg = 20.0, tilt_deg = 7.0;
std::vector<uint16_t> image(nx * ny, 5);
WriteThirdPartyDataFile("two_theta_000001.h5", image, 2, ny, nx);
{
HDF5File master("two_theta_master.h5");
HDF5Group entry(master, "entry");
entry.SaveScalar("definition", "NXmx");
HDF5Group instrument(entry, "instrument");
HDF5Group beam(instrument, "beam");
beam.SaveScalar("incident_wavelength", 0.6889)->Units("angstrom");
HDF5Group transformations(instrument, "transformations");
// Outermost first in the file, innermost first along the chain: det_z -> two_theta -> tilt
transformations.SaveVector("tilt", std::vector<double>{tilt_deg})
->Transformation("deg", ".", "detector", "", "rotation", {0, 1, 0});
transformations.SaveVector("two_theta", std::vector<double>{two_theta_deg})
->Transformation("deg", "/entry/instrument/transformations/tilt",
"detector", "", "rotation", {-1, 0, 0});
transformations.SaveVector("det_z", std::vector<double>{160.0})
->Transformation("mm", "/entry/instrument/transformations/two_theta",
"detector", "", "translation", {0, 0, 1});
HDF5Group detector(instrument, "detector");
detector.SaveScalar("depends_on", "/entry/instrument/transformations/det_z");
detector.SaveScalar("description", "PILATUS 2M");
detector.SaveScalar("beam_center_x", 4.0)->Units("pixels");
detector.SaveScalar("beam_center_y", 3.0)->Units("pixels");
detector.SaveScalar("distance", 0.160)->Units("m");
detector.SaveScalar("x_pixel_size", 0.172)->Units("mm");
detector.SaveScalar("y_pixel_size", 0.172)->Units("mm");
detector.SaveScalar("sensor_thickness", 0.32)->Units("mm");
detector.SaveScalar("count_time", 0.2);
detector.SaveScalar("saturation_value", static_cast<int64_t>(65535));
HDF5Group data(entry, "data");
data.ExternalLink("two_theta_000001.h5", "/data", "data_000001");
}
DiffractionGeometry geom;
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("two_theta_master.h5"));
geom = reader.GetDataset()->experiment.GetDiffractionGeometry();
}
// The chain as it stands in the internal frame: McStas (-1,0,0) is internal (1,0,0) and McStas
// (0,1,0) is internal (0,-1,0), and the outer rotation multiplies on the left.
const auto to_rad = [](double deg) { return static_cast<float>(deg * PI / 180.0); };
const RotMatrix expected = RotMatrix(to_rad(tilt_deg), {0, -1, 0})
* RotMatrix(to_rad(two_theta_deg), {1, 0, 0});
for (int64_t column = 0; column < 3; column++)
CHECK((geom.GetDetectorMatrix().Column(column) - expected.Column(column)).Length() < 1e-5f);
// Distance and beam centre are the ones the file states: the arm turns the detector about the
// sample and moves neither.
CHECK(geom.GetDetectorDistance_mm() == Catch::Approx(160.0));
CHECK(geom.GetBeamX_pxl() == Catch::Approx(4.0));
CHECK(geom.GetBeamY_pxl() == Catch::Approx(3.0));
// And the beam centre pixel is now that far from the beam - the whole point of a 2theta arm.
CHECK(geom.TwoTheta_rad(4.0f, 3.0f) * 180.0f / PI
== Catch::Approx(angle_deg(expected * Coord(0, 0, 1), Coord(0, 0, 1))));
CHECK(geom.TwoTheta_rad(4.0f, 3.0f) * 180.0f / PI > two_theta_deg);
remove("two_theta_000001.h5");
remove("two_theta_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// The same file with the arm parked at zero: the chain is there, nothing in it turns, and the
// geometry must be exactly the square-on one. This is nearly every file, so it has to cost nothing.
TEST_CASE("JFJochReader_DetectorTwoThetaZeroIsSquareOn", "[HDF5][Full]") {
const hsize_t nx = 8, ny = 6;
std::vector<uint16_t> image(nx * ny, 5);
WriteThirdPartyDataFile("two_theta_zero_000001.h5", image, 2, ny, nx);
{
HDF5File master("two_theta_zero_master.h5");
HDF5Group entry(master, "entry");
entry.SaveScalar("definition", "NXmx");
HDF5Group instrument(entry, "instrument");
HDF5Group beam(instrument, "beam");
beam.SaveScalar("incident_wavelength", 0.6889)->Units("angstrom");
HDF5Group transformations(instrument, "transformations");
transformations.SaveVector("two_theta", std::vector<double>{0.0})
->Transformation("deg", ".", "detector", "", "rotation", {-1, 0, 0});
transformations.SaveVector("det_z", std::vector<double>{160.0})
->Transformation("mm", "/entry/instrument/transformations/two_theta",
"detector", "", "translation", {0, 0, 1});
HDF5Group detector(instrument, "detector");
detector.SaveScalar("depends_on", "/entry/instrument/transformations/det_z");
detector.SaveScalar("description", "PILATUS 2M");
detector.SaveScalar("beam_center_x", 4.0)->Units("pixels");
detector.SaveScalar("beam_center_y", 3.0)->Units("pixels");
detector.SaveScalar("distance", 0.160)->Units("m");
detector.SaveScalar("x_pixel_size", 0.172)->Units("mm");
detector.SaveScalar("y_pixel_size", 0.172)->Units("mm");
detector.SaveScalar("sensor_thickness", 0.32)->Units("mm");
detector.SaveScalar("count_time", 0.2);
detector.SaveScalar("saturation_value", static_cast<int64_t>(65535));
HDF5Group data(entry, "data");
data.ExternalLink("two_theta_zero_000001.h5", "/data", "data_000001");
}
DiffractionGeometry geom;
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("two_theta_zero_master.h5"));
geom = reader.GetDataset()->experiment.GetDiffractionGeometry();
}
CHECK(geom.GetPoniRot1_rad() == 0.0f);
CHECK(geom.GetPoniRot2_rad() == 0.0f);
CHECK(geom.GetPoniRot3_rad() == 0.0f);
CHECK(geom.TwoTheta_rad(4.0f, 3.0f) == 0.0f);
remove("two_theta_zero_000001.h5");
remove("two_theta_zero_master.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// A file this system wrote states its own PONI angles twice: as the three scalars the reader takes
// them from, and as three rotations in the detector's depends_on chain. Following the chain must
// therefore skip them - applied on top of the scalars they would tilt the detector twice, which is
// how a correct 2theta reader breaks every tilted file this system has ever written. A test that only
// wrote an untilted detector could not see it.
TEST_CASE("JFJochReader_DetectorChainDoesNotDoubleTheTilt", "[HDF5][Full]") {
const float rot1 = 0.031f, rot2 = -0.047f, rot3 = 0.019f;
DiffractionExperiment x(DetJF(1));
x.ImagesPerTrigger(2).OverwriteExistingFiles(true).FilePrefix("test_ponichain");
x.BeamX_pxl(100).BeamY_pxl(200).DetectorDistance_mm(150)
.IncidentEnergy_keV(WVL_1A_IN_KEV).PixelSigned(false).BitDepthImage(16)
.FrameTime(std::chrono::microseconds(500), std::chrono::microseconds(10));
x.PoniRot1_rad(rot1).PoniRot2_rad(rot2).PoniRot3_rad(rot3);
RegisterHDF5Filter();
std::vector<uint16_t> image(x.GetPixelsNum(), 0);
StartMessage start_message;
x.FillMessage(start_message);
FileWriter file_set(start_message);
DataMessage message{};
for (int i = 0; i < x.GetImageNum(); i++) {
message.image = CompressedImage(image, x.GetXPixelsNum(), x.GetYPixelsNum());
message.number = i;
REQUIRE_NOTHROW(file_set.WriteHDF5(message));
}
EndMessage end_message;
end_message.max_image_number = x.GetImageNum();
file_set.WriteHDF5(end_message);
file_set.Finalize();
DiffractionGeometry geom;
{
JFJochHDF5Reader reader;
REQUIRE_NOTHROW(reader.ReadFile("test_ponichain_master.h5"));
geom = reader.GetDataset()->experiment.GetDiffractionGeometry();
}
CHECK(geom.GetPoniRot1_rad() == Catch::Approx(rot1).margin(1e-6));
CHECK(geom.GetPoniRot2_rad() == Catch::Approx(rot2).margin(1e-6));
CHECK(geom.GetPoniRot3_rad() == Catch::Approx(rot3).margin(1e-6));
for (int64_t column = 0; column < 3; column++)
CHECK((geom.GetDetectorMatrix().Column(column)
- PoniRotMatrix(rot1, rot2, rot3).Column(column)).Length() < 1e-5f);
remove("test_ponichain_master.h5");
remove("test_ponichain_data_000001.h5");
REQUIRE(H5Fget_obj_count(H5F_OBJ_ALL, H5F_OBJ_ALL) == 0);
}
// A miniCBF header states more about the instrument than the "# " lines do: the CBF template block
// some beamlines write carries a full imgCIF axis table, saying which laboratory direction the image's
// columns and rows run along and which the spindle turns about. The reader assumed all three, and two
// instruments in the corpus are not what it assumed - one stores its image a quarter turn round, the
// other turns about the VERTICAL. Either way the spindle ends up 90 degrees from the image, which is
// not a sign and so is beyond the axis-sign rescue; both indexed nothing at all.
namespace {
// Two frames of a sweep whose pixels are all zero. Every delta of a zero image is zero, so the
// byte-offset stream is one 0x00 per pixel - which is a valid stream and enough to open a sweep.
void WriteMiniCBFSweep(const std::string &prefix, const std::string &header_body,
int64_t nx, int64_t ny) {
for (int frame = 1; frame <= 2; frame++) {
std::ostringstream head;
head << "###CBF: VERSION 1.5\n_array_data.header_convention \"PILATUS_1.2\"\n"
<< "_array_data.header_contents\n"
<< "# Detector: PILATUS3 6M, S/N 60-0119\n"
<< "# Pixel_size 172e-6 m x 172e-6 m\n"
<< "# Silicon sensor, thickness 0.000450 m\n"
<< "# Exposure_time 0.1 s\n# Exposure_period 0.1 s\n# Count_cutoff 768595 counts\n"
<< "# Wavelength 0.96864 A\n# Detector_distance 0.33161 m\n"
<< "# Beam_xy (12.00, 8.00) pixels\n"
<< "# Start_angle " << (frame - 1) * 0.1 << " deg.\n# Angle_increment 0.1000 deg.\n"
<< "# Omega " << (frame - 1) * 0.1 << " deg.\n# Omega_increment 0.1000 deg.\n"
<< "# Phi 0.0000 deg.\n# Phi_increment 0.0000 deg.\n"
<< "# Chi 0.0000 deg.\n# Chi_increment 0.0000 deg.\n"
<< header_body
<< "_array_data.data\n--CIF-BINARY-FORMAT-SECTION--\n"
<< "Content-Type: application/octet-stream;\n"
<< " conversions=\"x-CBF_BYTE_OFFSET\"\n"
<< "Content-Transfer-Encoding: BINARY\n"
<< "X-Binary-Size: " << nx * ny << "\n"
<< "X-Binary-Element-Type: \"signed 32-bit integer\"\n"
<< "X-Binary-Number-of-Elements: " << nx * ny << "\n"
<< "X-Binary-Size-Fastest-Dimension: " << nx << "\n"
<< "X-Binary-Size-Second-Dimension: " << ny << "\n\n";
std::ostringstream name;
name << prefix << "_" << std::setfill('0') << std::setw(4) << frame << ".cbf";
std::ofstream f(name.str(), std::ios::binary);
const std::string text = head.str();
f.write(text.data(), static_cast<std::streamsize>(text.size()));
f.write(reinterpret_cast<const char *>(minicbf::BINARY_SEPARATOR),
sizeof(minicbf::BINARY_SEPARATOR));
const std::vector<char> zeros(static_cast<size_t>(nx * ny), 0);
f.write(zeros.data(), static_cast<std::streamsize>(zeros.size()));
}
}
void RemoveMiniCBFSweep(const std::string &prefix) {
for (int frame = 1; frame <= 2; frame++) {
std::ostringstream name;
name << prefix << "_" << std::setfill('0') << std::setw(4) << frame << ".cbf";
remove(name.str().c_str());
}
}
// The axis table in the form these headers write it, several tags to a line.
std::string AxisTable(const std::string &rows, int64_t nx, int64_t ny) {
return "loop_\n_axis.id\n_axis.type\n_axis.equipment\n_axis.depends_on\n"
"_axis.vector[1] _axis.vector[2] _axis.vector[3]\n"
"_axis.offset[1] _axis.offset[2] _axis.offset[3]\n"
+ rows +
"loop_\n_array_structure_list.array_id\n_array_structure_list.index\n"
"_array_structure_list.dimension\n_array_structure_list.precedence\n"
"_array_structure_list.direction\n_array_structure_list.axis_set_id\n"
"ARRAY1 1 " + std::to_string(nx) + " 1 increasing ELEMENT_X\n"
"ARRAY1 2 " + std::to_string(ny) + " 2 increasing ELEMENT_Y\n"
"loop_\n_array_structure_list_axis.axis_set_id\n_array_structure_list_axis.axis_id\n"
"_array_structure_list_axis.displacement\n_array_structure_list_axis.displacement_increment\n"
"ELEMENT_X ELEMENT_X 0.0 0.1720\nELEMENT_Y ELEMENT_Y 0.0 0.1720\n";
}
}
TEST_CASE("JFJochCBFReader_AxisTableStatesTheMounting", "[HDF5][Full]") {
const int64_t nx = 24, ny = 16;
// A header that states nothing: the assumption, and the behaviour of nearly every file there is.
SECTION("no table, no hint - the assumption stands") {
WriteMiniCBFSweep("cbfaxis_plain", "# Detector_2theta 0.0000 deg.\n# Oscillation_axis OMEGA\n",
nx, ny);
JFJochCBFReader reader;
REQUIRE_NOTHROW(reader.ReadFiles("cbfaxis_plain_0001.cbf"));
const auto x = reader.GetDataset()->experiment;
REQUIRE(x.GetGoniometer().has_value());
CHECK((x.GetGoniometer()->GetAxis() - Coord(-1, 0, 0)).Length() < 1e-6f);
CHECK(x.GetDetectorSetup().GetImageOrientation().IsIdentity());
CHECK(x.GetDiffractionGeometry().GetPoniRot2_rad() == 0.0f);
reader.Close();
RemoveMiniCBFSweep("cbfaxis_plain");
}
// A spindle that turns about the VERTICAL, with the image mounted the usual way round. imgCIF Y is
// up and the internal frame's y is down, so the stated (0,1,0) is internal (0,-1,0) - and NOT the
// (-1,0,0) that was assumed, which is 90 degrees away and indexes nothing.
SECTION("vertical spindle, standard image") {
WriteMiniCBFSweep("cbfaxis_vert",
"# Detector_2theta 0.0000 deg.\n# Oscillation_axis X.CW +SLOW\n"
+ AxisTable("GON_OMEGA rotation goniometer . 0 1 0 . . .\n"
"DET_Z translation detector . 0 0 -1 0 0 0\n"
"ELEMENT_X translation detector DET_Z 1 0 0 -1 1 0\n"
"ELEMENT_Y translation detector ELEMENT_X 0 -1 0 0 0 0\n", nx, ny),
nx, ny);
JFJochCBFReader reader;
REQUIRE_NOTHROW(reader.ReadFiles("cbfaxis_vert_0001.cbf"));
const auto x = reader.GetDataset()->experiment;
REQUIRE(x.GetGoniometer().has_value());
CHECK((x.GetGoniometer()->GetAxis() - Coord(0, -1, 0)).Length() < 1e-6f);
// The image itself is standard, so nothing about it is turned - the axis was the whole error.
CHECK(x.GetDetectorSetup().GetImageOrientation().IsIdentity());
reader.Close();
RemoveMiniCBFSweep("cbfaxis_vert");
}
// The same vertical spindle, stated only by the "+SLOW" token, which is all a header with no axis
// table says. Two datasets from that instrument are in this state.
SECTION("vertical spindle from the +SLOW token alone") {
WriteMiniCBFSweep("cbfaxis_slow", "# Detector_2theta 0.0000 deg.\n# Oscillation_axis X.CW +SLOW\n",
nx, ny);
JFJochCBFReader reader;
REQUIRE_NOTHROW(reader.ReadFiles("cbfaxis_slow_0001.cbf"));
const auto x = reader.GetDataset()->experiment;
REQUIRE(x.GetGoniometer().has_value());
CHECK((x.GetGoniometer()->GetAxis() - Coord(0, -1, 0)).Length() < 1e-6f);
reader.Close();
RemoveMiniCBFSweep("cbfaxis_slow");
}
// An image stored a quarter turn round, on a detector swung out to 30 degrees. The two are read
// together or not at all: the arm turns about a laboratory axis, and which way that runs across
// the stored image is exactly what the mounting says.
SECTION("quarter-turned image on a swung arm") {
WriteMiniCBFSweep("cbfaxis_turn",
"# Detector_2theta 30.0000 deg.\n# Oscillation_axis OMEGA\n"
+ AxisTable("GON_OMEGA rotation goniometer . 1 0 0 . . .\n"
"DET_2THETA rotation detector . 1 0 0 . . .\n"
"DET_Z translation detector DET_2THETA 0 0 -1 0 0 0\n"
"ELEMENT_X translation detector DET_Z 0 1 0 -1 1 0\n"
"ELEMENT_Y translation detector ELEMENT_X 1 0 0 0 0 0\n", nx, ny),
nx, ny);
JFJochCBFReader reader;
REQUIRE_NOTHROW(reader.ReadFiles("cbfaxis_turn_0001.cbf"));
const auto x = reader.GetDataset()->experiment;
REQUIRE(x.GetGoniometer().has_value());
CHECK((x.GetGoniometer()->GetAxis() - Coord(1, 0, 0)).Length() < 1e-6f);
// fast = imgCIF (0,1,0) = internal (0,-1,0), slow = imgCIF (1,0,0) = internal (1,0,0)
CHECK(x.GetDetectorSetup().GetImageOrientation() == DetectorOrientation(false, 3));
// and the arm turns about its own stated axis, internal +x, by the stated 30 degrees
const auto geom = x.GetDiffractionGeometry();
const RotMatrix expected = RotMatrix(static_cast<float>(30.0 * PI / 180.0), {1, 0, 0})
* DetectorOrientation(false, 3).Matrix();
for (int64_t column = 0; column < 3; column++)
CHECK((geom.GetDetectorMatrix().Column(column) - expected.Column(column)).Length() < 1e-5f);
reader.Close();
RemoveMiniCBFSweep("cbfaxis_turn");
}
}
Reference in New Issue View Git Blame Copy Permalink