diff --git a/common/DiffractionExperiment.cpp b/common/DiffractionExperiment.cpp index 8f51fbb7..e257d871 100644 --- a/common/DiffractionExperiment.cpp +++ b/common/DiffractionExperiment.cpp @@ -660,8 +660,10 @@ void DiffractionExperiment::FillMessage(StartMessage &message) const { message.image_size_x = GetXPixelsNum(); message.image_size_y = GetYPixelsNum(); message.mirror_y = IsDetectorMirroredY(); - message.saturation_value = GetSaturationLimit() - 1; - message.error_value = GetUnderflow(); + message.saturation_value = SaturationValueFromLimit(GetSaturationLimit()); + // The marker actually stored in the pixels: UINTx_MAX unsigned, INTx_MIN signed. GetUnderflow() + // was written here, which is -1 for an unsigned image and so matches no pixel it can contain. + message.error_value = GetImageFillValue(); message.frame_time = std::chrono::duration(GetImageTime()).count(); message.count_time = std::chrono::duration(GetImageCountTime()).count(); message.number_of_images = GetImageNum(); @@ -677,7 +679,13 @@ void DiffractionExperiment::FillMessage(StartMessage &message) const { // readout depth is a constant of the detector and of no use to a data consumer, so reporting it // here would only mislead: it differs from the image depth exactly when summation widens the // image, which is when getting it wrong does damage. - message.bit_depth_readout = message.bit_depth_image; + // Only for unsigned images. DIALS remaps 2^bit_depth_readout-1 and -2 to -1 and -2 whenever the + // field is present, without looking at the pixel type: for an unsigned image those land below + // underload_value and are masked, which is exactly what we want, but for a signed one they land + // INSIDE the trusted range and a clipped pixel is integrated as a count of -2. Signed images + // read correctly without it; unsigned 32-bit cannot be read at all without it. + if (!IsPixelSigned()) + message.bit_depth_readout = message.bit_depth_image; // Lowest valid value, so a reader can tell the error marker from data. Unsigned images use 0; // signed ones reserve INTx_MIN as the marker (GetUnderflow), so the lowest valid is one above. message.underload_value = IsPixelSigned() ? GetUnderflow() + 1 : 0; diff --git a/common/DiffractionExperiment.h b/common/DiffractionExperiment.h index 79b61341..84d35205 100644 --- a/common/DiffractionExperiment.h +++ b/common/DiffractionExperiment.h @@ -17,6 +17,17 @@ #include "../image_analysis/spot_finding/SpotFindingSettings.h" #include "DatasetSettings.h" #include "ROIMap.h" + +// Jungfraujoch's saturation limit is EXCLUSIVE: it is the first value that is no longer a real +// count. A saturated pixel is set to the type maximum, which IS that limit when no detector cutoff +// is configured, so the highest value that is still a measurement is one below it. +// NXmx saturation_value is the opposite convention - the highest value that IS valid - and the XDS +// OVERLOAD parameter and the DIALS trusted_range read it the same way, inclusively. So the two are +// exactly one apart. Convert with these both ways round: the pair had drifted apart, the write side +// subtracting one and the read side not adding it back, which cost a count on every +// write-read-write cycle and compounded without bound. +constexpr int64_t SaturationValueFromLimit(int64_t limit) { return limit - 1; } +constexpr int64_t SaturationLimitFromValue(int64_t saturation_value) { return saturation_value + 1; } #include "InstrumentMetadata.h" #include "ImageFormatSettings.h" #include "DetectorSettings.h" diff --git a/docs/CHANGELOG.md b/docs/CHANGELOG.md index 067b1e34..81cdf0c6 100644 --- a/docs/CHANGELOG.md +++ b/docs/CHANGELOG.md @@ -14,6 +14,9 @@ This is an UNSTABLE release. It includes many experimental features, as well as * rugnux: the detector geometry is also logged in XDS's convention (`ORGX`/`ORGY`, detector axis vectors, rotation axis), so it can be compared directly with an XDS refinement. * HDF5: a data file missing next to a VDS master now reads as the error-pixel marker instead of zero counts, so those frames are masked rather than silently integrated as blank. * HDF5: a still or grid scan recorded with a Smargon head position is no longer read back as a single image. +* HDF5 and image stream: `saturation_value` no longer loses a count each time a file is read and written again; a reprocessed `_process.h5` kept shrinking its own trusted range. +* HDF5: `error_value` now reports the marker an unsigned image actually stores (`UINTx_MAX`); it was written as `-1`, which no unsigned pixel can hold. +* HDF5: `bit_depth_readout` is written for unsigned images only. On a signed image it made DIALS read saturated pixels as a valid count of -2. * The writer refuses a stream whose start message declares a different pixel format than its images carry, instead of writing a master that does not describe its own data. * HDF5: `module_offset` is written as a float with a proper unit vector, and every transformation offset declares `offset_units`, so a reader does not fall back to the axis's own units - degrees on a rotation - when interpreting a length. * The image stream can carry the sample transformation chain (`transformations`, in the END message) in mounting order, so a goniometer axis, the Smargon chi/phi and a grid stage are described together and unambiguously. It is optional - a producer that does not send it gets the same chain built by the writer - so no metadata version change is needed. The rotation width and end angles are derived by the writer rather than sent. diff --git a/docs/HDF5.md b/docs/HDF5.md index f2922756..7e16bb8f 100644 --- a/docs/HDF5.md +++ b/docs/HDF5.md @@ -161,7 +161,7 @@ File-level HDF5 attributes `file_name`, `file_time`, `HDF5_Version` are also set | `x_pixel_size`, `y_pixel_size` | base | m | | `serial_number` | base | | | `bit_depth_readout` | NXmx | bit depth of the **stored image**, not of the detector electronics - see below | -| `saturation_value` | NXmx | highest valid value | +| `saturation_value` | NXmx | highest valid value. Read inclusively by NXmx, by the DIALS `trusted_range` and by the XDS `OVERLOAD` parameter; a saturated pixel carries the value one above it | | `underload_value` | NXmx | lowest valid value: `0` for an unsigned image, `INTx_MIN + 1` for a signed one | | `flatfield_applied` | NXmx | | | `pixel_mask`, `pixel_mask_applied` | NXmx | `pixel_mask` is `[y, x]`, hard-linked from `detectorSpecific/pixel_mask` | @@ -445,6 +445,20 @@ group for compatibility with existing tooling: | `detector_number` | | detector identifier (Dectris convention) | | `mirror_y` (in `detectorSpecific`) | | whether the stored image is mirrored in Y relative to the raw readout; true is the MX convention (row 0 at the top) | | `error_value` | | masked/error pixel sentinel: `UINTx_MAX` unsigned, `INTx_MIN` signed (NXmx has no equivalent; `underload_value` is written too, and is one above it) | + +NeXus has no concept of a fill or no-data value — it expects bad pixels to be flagged in `pixel_mask`, +which Jungfraujoch also writes. The in-band sentinel above is a DECTRIS compatibility convention: +SIMPLON specifies that masked pixels are flagged with `2^bit_depth_image - 1`. + +For an **unsigned** image the sentinel and the saturation code are the same value, so a saturated +pixel and a masked one cannot be told apart — the FPGA collapses both onto `UINTx_MAX`. Signed +images keep them separate: `INTx_MIN` is the marker, `INTx_MAX` is saturation. + +`bit_depth_readout` is written for **unsigned** images only. DIALS remaps the top two codes of +`2^bit_depth_readout` to `-1` and `-2` whenever the field is present, regardless of the pixel type: +for an unsigned image those land below `underload_value` and are correctly masked, but for a signed +one they land inside the trusted range and a saturated pixel would be integrated as a count of `-2`. +Signed images are read correctly without the field; unsigned 32-bit cannot be read at all without it. | `bit_depth_image` | | stored image bit depth; equal to `bit_depth_readout` (DECTRIS convention, not NXmx) | | `acquisition_type` | | always `triggered` (Dectris convention) | | `jungfrau_conversion_applied` | | JUNGFRAU photon/keV conversion applied | diff --git a/reader/HDF5MetadataSource.cpp b/reader/HDF5MetadataSource.cpp index d064078b..03eb1405 100644 --- a/reader/HDF5MetadataSource.cpp +++ b/reader/HDF5MetadataSource.cpp @@ -748,7 +748,8 @@ HDF5MetadataSource::OpenResult HDF5MetadataSource::Open(const std::string &filen detector.SensorThickness_um(master_file->GetFloat("/entry/instrument/detector/sensor_thickness") * 1e6); if (master_file->Exists("/entry/instrument/detector/sensor_material")) detector.SensorMaterial(master_file->GetString("/entry/instrument/detector/sensor_material")); - detector.SaturationLimit(master_file->GetInt("/entry/instrument/detector/saturation_value")); + detector.SaturationLimit(SaturationLimitFromValue( + master_file->GetInt("/entry/instrument/detector/saturation_value"))); // The reader hands every image out as signed int32 whatever the file stored (see PixelSigned // below), so that is the container depth the rest of the code has to see. DetectorSetup defaults // DECTRIS to 16 bits and GetByteDepthImage() prefers the detector's value over the image diff --git a/reader/JFJochReaderImage.cpp b/reader/JFJochReaderImage.cpp index 57428ca6..2d22b00a 100644 --- a/reader/JFJochReaderImage.cpp +++ b/reader/JFJochReaderImage.cpp @@ -226,7 +226,7 @@ void JFJochReaderImage::AddImage(const JFJochReaderImage &other) { if (sum < MIN_REAL_PXL_VALUE) [[unlikely]] { image[i] = ERROR_PXL_VALUE; error_pixel.emplace(static_cast(i)); - } else if (sum > dataset->experiment.GetSaturationLimit()) [[unlikely]] { + } else if (sum >= dataset->experiment.GetSaturationLimit()) [[unlikely]] { image[i] = SATURATED_PXL_VALUE; saturated_pixel.emplace(static_cast(i)); } else { diff --git a/receiver/JFJochReceiverLite.cpp b/receiver/JFJochReceiverLite.cpp index 06b8ab7a..d642ea47 100644 --- a/receiver/JFJochReceiverLite.cpp +++ b/receiver/JFJochReceiverLite.cpp @@ -204,7 +204,7 @@ void JFJochReceiverLite::Configure(const StartMessage &msg) { experiment.Detector().SerialNumber(msg.detector_serial_number); experiment.Detector().SensorMaterial(msg.sensor_material); experiment.Detector().SensorThickness_um(msg.sensor_thickness * 1e6); - experiment.Detector().SaturationLimit(msg.saturation_value); + experiment.Detector().SaturationLimit(SaturationLimitFromValue(msg.saturation_value)); // Images are forwarded byte-for-byte, so the stream's own image_dtype - not anything configured // locally - decides both the width and the sign the outgoing metadata must declare. Taking only // the width used to leave a signed stream declared unsigned in NXmx. diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index 8fc8c88a..90b44269 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -985,7 +985,11 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b // mapping and opens the source itself. const auto stored = reader_.GetStoredPixelFormat(); start_message.bit_depth_image = stored.bit_depth; - start_message.bit_depth_readout = stored.bit_depth; + // Unsigned only, for the reason FillMessage gives: DIALS remaps the top two codes whenever the + // field is present, and on a signed image those land inside the trusted range. + start_message.bit_depth_readout.reset(); + if (!stored.is_signed) + start_message.bit_depth_readout = stored.bit_depth; start_message.pixel_signed = stored.is_signed; start_message.error_value = stored.is_signed ? -(int64_t(1) << (stored.bit_depth - 1)) // INTx_MIN diff --git a/tests/DiffractionExperimentTest.cpp b/tests/DiffractionExperimentTest.cpp index 120b502e..24f32900 100644 --- a/tests/DiffractionExperimentTest.cpp +++ b/tests/DiffractionExperimentTest.cpp @@ -1133,3 +1133,49 @@ TEST_CASE("DiffractionExperiment_ImagesPerFile_Defaults", "[DiffractionExperimen CHECK(x.GetImagesPerFile() == 5000); } } + +// The saturation limit is exclusive inside Jungfraujoch and inclusive in the file, so the two sides +// must convert. They stopped agreeing once, and the value then lost a count on every +// write-read-write cycle; these pin the pairing and what each stored type declares. +TEST_CASE("DiffractionExperiment_SaturationAndErrorMarkers", "[DiffractionExperiment]") { + SECTION("the conversions are inverses") { + for (int64_t limit : {256L, 65536L, 32768L, 2147483648L}) { + CHECK(SaturationLimitFromValue(SaturationValueFromLimit(limit)) == limit); + CHECK(SaturationValueFromLimit(limit) == limit - 1); + } + } + + struct Case { int64_t bits; bool is_signed; int64_t saturation; int64_t underload; int64_t error; }; + // saturation_value is the highest value that is still a count; error is the marker actually + // stored in the pixels. Unsigned images reserve the top code, signed ones the bottom. + const std::vector cases = { + {16, false, UINT16_MAX - 1, 0, UINT16_MAX}, + {32, false, UINT32_MAX - 1, 0, UINT32_MAX}, + {16, true, INT16_MAX - 1, INT16_MIN + 1, INT16_MIN}, + {32, true, INT32_MAX - 1, INT32_MIN + 1, INT32_MIN}, + }; + + for (const auto &c: cases) { + DiffractionExperiment x(DetJF(1)); + x.BitDepthImage(c.bits).PixelSigned(c.is_signed); + StartMessage message; + x.FillMessage(message); + + CAPTURE(c.bits, c.is_signed); + CHECK(message.saturation_value == c.saturation); + CHECK(message.saturation_value == SaturationValueFromLimit(x.GetSaturationLimit())); + REQUIRE(message.underload_value.has_value()); + CHECK(message.underload_value.value() == c.underload); + // The marker the pixels really carry - GetUnderflow() was written here, and it is -1 for an + // unsigned image, which matches no pixel such an image can hold. + REQUIRE(message.error_value.has_value()); + CHECK(message.error_value.value() == c.error); + CHECK(message.error_value.value() == x.GetImageFillValue()); + + // DIALS remaps the top two codes whenever bit_depth_readout is present, without looking at + // the pixel type. For an unsigned image they fall below underload_value and are masked, as + // intended; for a signed one they land inside the trusted range and a saturated pixel would + // be integrated as a count of -2. Unsigned 32-bit cannot be read at all without the field. + CHECK(message.bit_depth_readout.has_value() == !c.is_signed); + } +} diff --git a/tests/JFJochReaderTest.cpp b/tests/JFJochReaderTest.cpp index 1c8393da..ac6f2369 100644 --- a/tests/JFJochReaderTest.cpp +++ b/tests/JFJochReaderTest.cpp @@ -3010,6 +3010,68 @@ TEST_CASE("JFJochReader_TransformationChain_SentAndBuilt", "[HDF5][Full]") { 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 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 diff --git a/viewer/JFJochHttpReader.cpp b/viewer/JFJochHttpReader.cpp index 95210bf2..f0111062 100644 --- a/viewer/JFJochHttpReader.cpp +++ b/viewer/JFJochHttpReader.cpp @@ -232,7 +232,7 @@ std::shared_ptr JFJochHttpReader::UpdateDataset_i() { DetectorSetup detector = DetDECTRIS(msg->start_message->image_size_x, msg->start_message->image_size_y, msg->start_message->detector_description, {}); detector.PixelSize_um(msg->start_message->pixel_size_x * 1e6); - detector.SaturationLimit(msg->start_message->saturation_value); + detector.SaturationLimit(SaturationLimitFromValue(msg->start_message->saturation_value)); detector.MinFrameTime(std::chrono::microseconds(0)); detector.MinCountTime(std::chrono::microseconds(0)); detector.ReadOutTime(std::chrono::microseconds (0));