diff --git a/image_analysis/bragg_prediction/BraggPrediction.cpp b/image_analysis/bragg_prediction/BraggPrediction.cpp index ac483c60..5e4460f3 100644 --- a/image_analysis/bragg_prediction/BraggPrediction.cpp +++ b/image_analysis/bragg_prediction/BraggPrediction.cpp @@ -4,9 +4,21 @@ #include #include "../../common/JFJochMath.h" +#include "../../common/Logger.h" #include "BraggPrediction.h" #include "../bragg_integration/SystematicAbsence.h" +void BraggPrediction::ReportOverflow(int predicted) { + if (overflow_reported) + return; + overflow_reported = true; + Logger("BraggPrediction").Warning( + "A frame predicted {} reflections but the buffer holds {} - the ones kept are whichever the GPU " + "wrote first, so this dataset does NOT process reproducibly. Results will differ between runs. " + "The cell is large enough that the prediction cap no longer fits it.", + predicted, kPredictionCapacity); +} + int BraggPrediction::TruncateToOutput(int count) { if (count <= kPredictionOutput) return count; diff --git a/image_analysis/bragg_prediction/BraggPrediction.h b/image_analysis/bragg_prediction/BraggPrediction.h index 60400b00..479dcf63 100644 --- a/image_analysis/bragg_prediction/BraggPrediction.h +++ b/image_analysis/bragg_prediction/BraggPrediction.h @@ -30,14 +30,23 @@ protected: const int max_reflections; std::vector reflections; + // A frame that predicts more than the buffer holds keeps an ARBITRARY subset of them: the GPU + // kernels claim slots with an atomicAdd, so which ones survive depends on block scheduling and + // changes from run to run. Say so, once per predictor, rather than let it pass silently - it is + // not a small effect (measured: a 3% run-to-run spread in the number of merged reflections, and + // every frame of that dataset overflowed). + void ReportOverflow(int predicted); + bool overflow_reported = false; + // Deterministically cap Calc's output at kPredictionOutput: if more were predicted, keep the ones // closest to the Ewald sphere (smallest excitation error), ties broken by hkl. Returns the kept // count. Below the cap it is a no-op. Call at the end of every Calc override. int TruncateToOutput(int count); public: // The prediction buffer holds up to kPredictionCapacity reflections so a strong lattice does not - // overflow it (the GPU kernels then fill it in a non-deterministic atomic order). Calc returns at - // most kPredictionOutput, the number that flows downstream and is serialized - kept low so the + // overflow it (the GPU kernels then fill it in a non-deterministic atomic order - ReportOverflow + // warns when that happens; a large unit cell, ~2.8e6 A^3, reaches ~40000 per frame and overflows on + // every one). Calc returns at most kPredictionOutput, the number that flows downstream and is serialized - kept low so the // per-image reflection list stays within the frame transport headroom. static constexpr int kPredictionCapacity = 20000; static constexpr int kPredictionOutput = 10000; diff --git a/image_analysis/bragg_prediction/BraggPredictionGPU.cu b/image_analysis/bragg_prediction/BraggPredictionGPU.cu index 6a606182..43ebefbf 100644 --- a/image_analysis/bragg_prediction/BraggPredictionGPU.cu +++ b/image_analysis/bragg_prediction/BraggPredictionGPU.cu @@ -148,9 +148,9 @@ namespace { int l = li - max_hkl; Reflection r{}; if (!compute_reflection(*kc, h, k, l, r)) return; - int pos = atomicAdd(counter, 1); + // See the rotation kernel: clamping the counter hides an overflow from the host. + const int pos = atomicAdd(counter, 1); if (pos < max_reflections) out[pos] = r; - else atomicSub(counter, 1); } inline KernelConsts BuildKernelConsts(const DiffractionExperiment &experiment, @@ -212,7 +212,10 @@ int BraggPredictionGPU::Calc(const DiffractionExperiment &experiment, cudaStreamSynchronize(stream); int count = *h_count.get(); - if (count > max_reflections) count = max_reflections; + if (count > max_reflections) { + ReportOverflow(count); + count = max_reflections; + } if (count == 0) return {}; diff --git a/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu b/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu index 5d317abe..c30d7b5f 100644 --- a/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu +++ b/image_analysis/bragg_prediction/BraggPredictionRotGPU.cu @@ -189,11 +189,11 @@ namespace { int n = compute_reflections_rot(*kc, h, k, l, r); for (int i = 0; i < n; ++i) { - int pos = atomicAdd(counter, 1); + // Do NOT clamp the counter back down on overflow: it then saturates at the capacity and the + // host cannot tell a full buffer from an overflowing one. Let it count the true total. + const int pos = atomicAdd(counter, 1); if (pos < max_reflections) out[pos] = r[i]; - else - atomicSub(counter, 1); } } @@ -289,7 +289,10 @@ int BraggPredictionRotGPU::Calc(const DiffractionExperiment &experiment, cudaStreamSynchronize(stream); int count = *h_count.get(); - if (count > max_reflections) count = max_reflections; + if (count > max_reflections) { + ReportOverflow(count); + count = max_reflections; + } if (count == 0) return 0; cudaMemcpyAsync(reflections.data(), d_out, sizeof(Reflection) * count, cudaMemcpyDeviceToHost, stream);