Take the serial passes out of the combine phase

Four passes over all ten million partials and one over all four million fulls,
each of them a single thread walking a fat struct: the snapshot of corr taken
before the pass filters, the zeta filter, the frame rejection, the extraction of
corr for the upload, and the scatter of the downloaded fulls arrays back into
Obs. All are independent per element, so all are chunked now, and the split is
gated on the amount of data rather than the core count.

Measured on the heaviest crystal, summed over the run: corr snapshot 0.14 s ->
0.05, zeta filter 0.20 -> 0.04, corr upload 0.16 -> 0.07, fulls scatter 0.27 ->
0.06.

Two things were tried in the same phase and are NOT here, both measured on the
way past. Spreading the combine over the other three GPUs: the whole
gpu_->Combine() call is 1.27 s over the entire run, so four cards could save
about a second at best, but the partials would have to live on every device -
1.7 GB of transfer at the 1.5 GB/s this machine gets - and the fulls would have
to be gathered back for the merge, which reads them where they are. The
replication alone costs more than the best case saves. And staging the fulls
download through a page-locked block, which reads like the textbook fix for
copies running at 1.5 GB/s, measured 0.76 s -> 0.80 s: whatever limits them
here, it is not the driver's bounce buffer.

Battery 9m23s, space group 21/24, no failures; the crystals that moved are the
two already known to sit on knife edges.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
This commit is contained in:
jungfrau
2026-08-16 05:59:16 -04:00
co-authored by Claude Opus 5
parent a1325637d2
commit 8ea3076f96
@@ -3182,19 +3182,29 @@ RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search) {
const bool pass_filters = (for_search && search_min_zeta > 0.0) || min_cc_for_image > 0.0;
if (pass_filters) {
corr_before_pass_filters.resize(partials.size());
for (size_t i = 0; i < partials.size(); ++i)
corr_before_pass_filters[i] = partials[i].corr;
ParallelChunks(static_cast<int>(partials.size()),
ThreadsForWork(partials.size(), nthreads), [&](int lo, int hi) {
for (int i = lo; i < hi; ++i) corr_before_pass_filters[i] = partials[i].corr;
});
}
bool corr_filtered = false;
// --- 2a. On the de-novo search pass only, drop observations whose Lorentz geometry is poor. ---
if (for_search && search_min_zeta > 0.0) {
int64_t n_dropped = 0;
for (auto &o : partials)
if (!(std::isfinite(o.zeta) && o.zeta >= search_min_zeta)) {
if (std::isfinite(o.corr) && o.corr > 0.0f) ++n_dropped;
o.corr = 0.0f;
std::atomic<int64_t> dropped = 0;
ParallelChunks(static_cast<int>(partials.size()),
ThreadsForWork(partials.size(), nthreads), [&](int lo, int hi) {
int64_t local = 0;
for (int i = lo; i < hi; ++i) {
Obs &o = partials[i];
if (!(std::isfinite(o.zeta) && o.zeta >= search_min_zeta)) {
if (std::isfinite(o.corr) && o.corr > 0.0f) ++local;
o.corr = 0.0f;
}
}
dropped += local;
});
const int64_t n_dropped = dropped.load();
corr_filtered = corr_filtered || n_dropped > 0;
if (n_dropped > 0)
logger.Info("Space-group search: ignoring {} observations with |zeta| < {:.2f} "
@@ -3212,8 +3222,11 @@ RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search) {
++n_rejected;
}
if (n_rejected > 0) {
for (auto &o : partials)
if (reject[o.frame]) o.corr = 0.0f;
ParallelChunks(static_cast<int>(partials.size()),
ThreadsForWork(partials.size(), nthreads), [&](int lo, int hi) {
for (int i = lo; i < hi; ++i)
if (reject[partials[i].frame]) partials[i].corr = 0.0f;
});
corr_filtered = true;
logger.Info("Rejected {} of {} frames correlating below {:.2f} with the merged reference",
n_rejected, n_frames, min_cc_for_image);
@@ -3226,7 +3239,10 @@ RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search) {
// left it doing nothing whenever the CPU combine was in use (the diagnostic dump).
if (gpu_active_ && corr_filtered) {
std::vector<float> corr(partials.size());
for (size_t i = 0; i < partials.size(); ++i) corr[i] = partials[i].corr;
ParallelChunks(static_cast<int>(partials.size()),
ThreadsForWork(partials.size(), nthreads), [&](int lo, int hi) {
for (int i = lo; i < hi; ++i) corr[i] = partials[i].corr;
});
gpu_->SetCorr(corr.data());
}
#endif
@@ -3280,7 +3296,8 @@ RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search) {
gpu_->GetFullsPxPy(fpx.data(), fpy.data());
gpu_->GetFullsVariance(fvb.data(), fvi.data());
if (scaled_fulls_on_gpu) gpu_->GetFullsCorr(fcorr.data());
for (int i = 0; i < nf; ++i) {
ParallelChunks(nf, ThreadsForWork(static_cast<size_t>(nf), nthreads), [&](int lo, int hi) {
for (int i = lo; i < hi; ++i) {
Obs &o = fulls[i];
o.h = fh[i]; o.k = fk[i]; o.l = fl[i];
o.I = fI[i]; o.sigma = fsig[i]; o.d = fd[i];
@@ -3289,6 +3306,7 @@ RotationScaleMerge::Result RotationScaleMerge::Run(bool for_search) {
o.image_number = fimg[i]; o.frame = fframe[i]; o.px = fpx[i]; o.py = fpy[i];
o.on_ice = fon[i]; o.group = fgroup[i];
}
});
logger.Info("3D combine{} (GPU): {} fulls", scaled_fulls_on_gpu ? " + scale-fulls" : "", nf);
combined_on_gpu = true;
}