diff --git a/rugnux/RugnuxDefaults.cpp b/rugnux/RugnuxDefaults.cpp index 8ef13bbb4..419228539 100644 --- a/rugnux/RugnuxDefaults.cpp +++ b/rugnux/RugnuxDefaults.cpp @@ -5,6 +5,9 @@ #include #include +#include + +#include "../common/ParallelFor.h" namespace { // Smoothing range for the per-frame scale G on rotation data (XDS DELPHI-like), in degrees. @@ -125,18 +128,20 @@ bool RugnuxDetectIceRings(std::optional file_detect_ice_rings, bool rotati return file_detect_ice_rings.value_or(rotation); } -RugnuxMeasuredGain RugnuxMeasureGain(JFJochReader &reader) { +RugnuxMeasuredGain RugnuxMeasureGain(JFJochReader &reader, size_t nthreads) { const auto dataset = reader.GetDataset(); const size_t n = reader.GetNumberOfImages(); const size_t frames = std::min(n, background_gain::GAIN_MAX_FRAMES); - std::vector levels; - std::vector histogram(background_gain::PEAK_HISTOGRAM_BINS, 0); - JFJochReaderRawImage raw; - std::vector buffer; - std::vector pixels; - for (size_t i = 0; i < frames; i++) { + // The frames are read and measured side by side; each keeps its own levels and histogram, which are + // put together in frame order below. + std::vector> frame_levels(frames); + std::vector> frame_histogram(frames); + ParallelFor(static_cast(frames), std::min(frames, nthreads), [&](int i) { + JFJochReaderRawImage raw; + std::vector buffer; + std::vector pixels; if (!reader.ReadRawImage(static_cast(i * n / frames), raw)) - continue; + return; const size_t w = raw.image.GetWidth(), h = raw.image.GetHeight(); pixels.resize(w * h); const uint8_t *ptr = raw.image.GetUncompressedPtr(buffer); @@ -147,11 +152,20 @@ RugnuxMeasuredGain RugnuxMeasureGain(JFJochReader &reader) { else if (raw.image.GetMode() == CompressedImageMode::Uint16) std::copy_n(reinterpret_cast(ptr), w * h, pixels.begin()); else - continue; - levels.push_back(background_gain::MeasureFrame(pixels, w, h, dataset->pixel_mask->GetMask(), - dataset->experiment.GetSaturationLimit())); + return; + frame_levels[i] = background_gain::MeasureFrame(pixels, w, h, dataset->pixel_mask->GetMask(), + dataset->experiment.GetSaturationLimit()); + frame_histogram[i].assign(background_gain::PEAK_HISTOGRAM_BINS, 0); background_gain::AddToHistogram(pixels, dataset->pixel_mask->GetMask(), - dataset->experiment.GetSaturationLimit(), histogram); + dataset->experiment.GetSaturationLimit(), frame_histogram[i]); + }); + std::vector levels; + std::vector histogram(background_gain::PEAK_HISTOGRAM_BINS, 0); + for (size_t i = 0; i < frames; i++) { + if (frame_levels[i]) + levels.push_back(*frame_levels[i]); + for (size_t b = 0; b < frame_histogram[i].size(); b++) + histogram[b] += frame_histogram[i][b]; } return {background_gain::PhotonPeak(histogram), background_gain::Gain(levels)}; } diff --git a/rugnux/RugnuxDefaults.h b/rugnux/RugnuxDefaults.h index 8d53f694b..9c8771b0b 100644 --- a/rugnux/RugnuxDefaults.h +++ b/rugnux/RugnuxDefaults.h @@ -84,7 +84,7 @@ struct RugnuxMeasuredGain { std::optional peaks; std::optional variance; }; -RugnuxMeasuredGain RugnuxMeasureGain(JFJochReader &reader); +RugnuxMeasuredGain RugnuxMeasureGain(JFJochReader &reader, size_t nthreads); // The gain the frames are analysed at, into everything that takes a photon count from them: the // integration (BraggIntegrationSettings::pixel_gain, which the defective-pixel tests read too) and the diff --git a/rugnux/rugnux_cli.cpp b/rugnux/rugnux_cli.cpp index 4c6cc6071..90bb92fca 100644 --- a/rugnux/rugnux_cli.cpp +++ b/rugnux/rugnux_cli.cpp @@ -2156,7 +2156,7 @@ int RugnuxInputRun::PrepareRugnuxInput() { const float unit = float_frames ? 1.0f / 16.0f : 1.0f; if (float_frames) xfel_reader.ApplyGain(unit); - const RugnuxMeasuredGain g = RugnuxMeasureGain(reader); + const RugnuxMeasuredGain g = RugnuxMeasureGain(reader, nthreads); if (float_frames) xfel_reader.ClearGain(); if (g.peaks) { diff --git a/viewer/JFJochProcessController.cpp b/viewer/JFJochProcessController.cpp index 3c3ca57a3..5ba76cbfb 100644 --- a/viewer/JFJochProcessController.cpp +++ b/viewer/JFJochProcessController.cpp @@ -83,7 +83,7 @@ void JFJochProcessController::run_(QString file_path, DiffractionExperiment expe JFJochReader &reader = *reader_ptr; // A CCD's pixel values are ADU: analysed at the gain its background shows, as rugnux does. if (reader_ptr == &marccd_reader || reader_ptr == &smv_reader) { - const RugnuxMeasuredGain g = RugnuxMeasureGain(reader); + const RugnuxMeasuredGain g = RugnuxMeasureGain(reader, static_cast(std::max(1, config.nthreads))); const std::optional gain = g.peaks ? std::optional(g.peaks->gain) : g.variance; ApplyRugnuxPixelGain(experiment, config.spot_finding, gain.value_or(1.0f)); }