Files
Jungfraujoch/reader/JFJochReader.cpp
T
leonarski_fandClaude Opus 5.5 aa3fa6b9c7 Reader: stream the frame data into the page cache ahead of the image loops
A cold run on a spinning disk waited on the disk twice over. The CBF header
scan read 256 kB from every frame on eight threads that each strode through
their own share of the sweep, so they drifted apart and the scan became a
seek storm (34 s for 2400 frames here); and after it, the pre-scan and the
first-pass indexing touch a few hundred frames and leave the disk idle until
the first image loop reads everything at seek-bound rates.

- ReadAhead (reader/): once the dataset is open, rugnux starts eight threads
  that read the data files - HDF5 data files (legacy, VDS or the integrated
  master) or the per-frame CBF/marCCD/SMV files - in 4 MB pieces taken
  strictly in order, into a throwaway buffer. One stream reads this disk at
  125 MB/s, eight in-order streams at 190 MB/s, 32 at 157 MB/s. It never gets
  more than a quarter of MemAvailable (GlobalMemoryStatusEx on Windows, 4 GiB
  where there is no figure) ahead of what ReadRawImage has handed out, so a
  dataset bigger than the cache does not evict its own start, and it stops
  with the reader. Plain ifstream reads: portable, no POSIX calls.
- Header scans (CBF, marCCD, SMV) hand the files out in order from an atomic
  counter (sweep::ForEachInOrder) instead of striding: 18 s -> 12 s for 2400
  cold CBF headers. The CBF header is first read with a 16 kB probe and again
  with the old 256 kB one only when the separator is not in it, so the parsed
  header is exactly what it was: 12 s -> 6 s.

Output unchanged: p.hkl, p.mtz and p_unmerged.mtz md5-identical to the
rc173 baseline on 6toc (CBF, 2400 frames, 6.0 GB) and 9q41 (HDF5 VDS, 900
frames, 5.1 GB), and on 6z9g (HDF5, 12.8 GB) to the unmodified branch; myob
(p.hkl p.mtz p_P1.mtz p_unmerged.mtz) md5-identical to the reference.

Measured cold (files evicted with POSIX_FADV_DONTNEED before every run),
same code without this commit vs with it, on a shared box (load 20-70, other
agents reading the same disk, so single runs scatter by +-20 s):
  6toc  wall 61.7/62.7 -> 49.3/49.8 s (clean pairs); all data resident
        after 62/51/50 -> 45/41/42 s
  9q41  wall 67.3 -> 57.8 s (clean pair); resident after 58/46/43 -> 48/37/35 s
  6z9g  resident after 81 -> 69 s
The first image loop can look slower with this in CBF runs: the old 256 kB
header probes pulled ~70% of the data in as kernel readahead, so the old
loop started warm - after a 34 s header scan instead of 10 s.
Warm (myob, NVMe, cached): 19.35/20.00 s without, 19.76-20.16 s with; the
read-ahead then only copies 9.3 GB out of the page cache, 0.44 s wall and
3.4 CPU-s measured standalone.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01D1G8gJVAy6gp1K5Dz3NE5C
2026-09-27 10:32:58 +02:00

148 lines
5.9 KiB
C++

// SPDX-FileCopyrightText: 2025 Filip Leonarski, Paul Scherrer Institute <filip.leonarski@psi.ch>
// SPDX-License-Identifier: GPL-3.0-only
#include "JFJochReader.h"
#include <future>
JFJochReader &JFJochReader::Experiment(const DiffractionExperiment &experiment) {
std::unique_lock ul(m);
default_experiment = experiment;
return *this;
}
void JFJochReader::SummationThread(int64_t image0, int64_t n_image, int64_t image_jump, JFJochReaderImage &image) {
std::vector<uint8_t> buffer;
DataMessage msg;
for (int64_t i = image0; i < n_image; i += image_jump) {
bool ret = LoadImage_i(dataset, msg, buffer, i, false);
if (ret) {
auto image_sum = std::make_shared<JFJochReaderImage>(msg, dataset);
{
std::unique_lock ul(summation_mutex);
image.AddImage(*image_sum);
}
}
}
}
std::shared_ptr<JFJochReaderImage> JFJochReader::LoadImage(int64_t image_number, int64_t summation_factor) {
// It would be a mess to load two images at the same time
// so loading is protected via mutex
// yet copying share_ptr pointer is atomic and needs no mutex protection
std::unique_lock ul(m);
std::vector<uint8_t> buffer;
DataMessage msg;
if (LoadImage_i(dataset, msg, buffer, image_number, true)) {
auto image = std::make_shared<JFJochReaderImage>(msg, dataset);
if (summation_factor > 4) {
int64_t nthread = std::min<int64_t>(summation_factor - 1, 8);
std::vector<std::future<void>> futures;
for (int i = 0; i < nthread; i++)
futures.emplace_back(std::async(std::launch::async,
&JFJochReader::SummationThread, this,
image_number + 1 + i,
image_number + summation_factor,
nthread,
std::ref(*image)));
for (auto &f: futures)
f.get();
} else if (summation_factor > 1) {
SummationThread(image_number + 1, image_number + summation_factor, 1, *image);
}
return image;
}
return {};
}
std::shared_ptr<JFJochReaderRawImage> JFJochReader::GetRawImage(int64_t image_number) {
auto ret = std::make_shared<JFJochReaderRawImage>();
if (!ReadRawImage(image_number, *ret))
return {};
return ret;
}
void JFJochReader::StartReadAhead() {
read_ahead_.reset();
const auto files = DataFiles();
if (!files.empty() && GetNumberOfImages() > 0)
read_ahead_ = std::make_unique<ReadAhead>(files, GetNumberOfImages());
}
void JFJochReader::SetStartMessage(const std::shared_ptr<JFJochReaderDataset> &val) {
std::unique_lock ul(m);
dataset = val;
}
std::shared_ptr<const JFJochReaderDataset> JFJochReader::GetDataset() const {
std::unique_lock ul(m);
return dataset;
}
void JFJochReader::UpdateGeomMetadata(const DiffractionExperiment &experiment) {
std::unique_lock ul(m);
if (!dataset)
return;
auto new_dataset = std::make_shared<JFJochReaderDataset>(*dataset);
// At the moment subset of options is limited to safe ones...need to change it in the future
new_dataset->experiment.BeamX_pxl(experiment.GetBeamX_pxl());
new_dataset->experiment.BeamY_pxl(experiment.GetBeamY_pxl());
new_dataset->experiment.DetectorDistance_mm(experiment.GetDetectorDistance_mm());
new_dataset->experiment.IncidentEnergy_keV(experiment.GetIncidentEnergy_keV());
new_dataset->experiment.PoniRot1_rad(experiment.GetDatasetSettings().GetPoniRot1_rad());
new_dataset->experiment.PoniRot2_rad(experiment.GetDatasetSettings().GetPoniRot2_rad());
new_dataset->experiment.PoniRot3_rad(experiment.GetDatasetSettings().GetPoniRot3_rad());
new_dataset->experiment.SetUnitCell(experiment.GetUnitCell());
new_dataset->experiment.SetSpaceGroup(experiment.GetGemmiSpaceGroup());
new_dataset->experiment.PolarizationFactor(experiment.GetPolarizationFactor());
new_dataset->experiment.Goniometer(experiment.GetGoniometer());
new_dataset->experiment.GridScan(experiment.GetGridScan());
new_dataset->experiment.ImportIndexingSettings(experiment.GetIndexingSettings());
new_dataset->experiment.ImportBraggIntegrationSettings(experiment.GetBraggIntegrationSettings());
new_dataset->experiment.DetectIceRings(experiment.IsDetectIceRings());
dataset = new_dataset;
}
void JFJochReader::UpdateUserMask(const std::vector<uint32_t> &mask) {
std::unique_lock ul(m);
if (!dataset)
return;
auto new_dataset = std::make_shared<JFJochReaderDataset>(*dataset);
// Copy-on-write: the mask is shared with the old snapshot, so edit a fresh copy, not in place.
auto new_mask = std::make_shared<PixelMask>(*dataset->pixel_mask);
new_mask->LoadUserMask(dataset->experiment, mask);
new_dataset->pixel_mask = new_mask;
dataset = new_dataset;
}
std::shared_ptr<JFJochReaderSpots> JFJochReader::ReadAllSpots(int64_t start_image, int64_t end_image,
int64_t stride) const {
if (start_image < 0)
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Start image must be non-negative");
if (start_image > end_image)
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Start image number is greater than end image number");
if (stride == 0)
throw JFJochException(JFJochExceptionCategory::InputParameterInvalid,
"Stride cannot be zero");
size_t nelems = (end_image - start_image) / stride + 1;
auto ret = std::make_shared<JFJochReaderSpots>();
ret->start_image = static_cast<int64_t>(start_image);
ret->stride = static_cast<int64_t>(stride);
ret->spots.reserve(nelems);
for (int i = 0; i < nelems; i++)
ret->spots.emplace_back(ReadSpots(start_image + i * stride));
return ret;
}