From 92877e5eb786e877380cfd462f2ca5895540bcb8 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Sat, 26 Sep 2026 22:45:48 +0200 Subject: [PATCH] Rugnux: build the unmerged MTZ beside the P1 cross-check merge The unmerged file reads the integration outcomes and the determined group, neither of which the P1 merge changes, except each image's mosaicity, which the batch headers carry and the merge rewrites. UnmergedMtz builds the file without it on a second thread, and SetUnmergedMtzMosaicity fills it in after the merge, so the file is the same bytes as before. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01D1G8gJVAy6gp1K5Dz3NE5C --- image_analysis/WriteReflections.cpp | 38 ++++++++++++++++++++++------- image_analysis/WriteReflections.h | 12 +++++++++ rugnux/Rugnux.cpp | 23 +++++++++++++++-- 3 files changed, 62 insertions(+), 11 deletions(-) diff --git a/image_analysis/WriteReflections.cpp b/image_analysis/WriteReflections.cpp index dff8b25a5..bd654a999 100644 --- a/image_analysis/WriteReflections.cpp +++ b/image_analysis/WriteReflections.cpp @@ -736,12 +736,11 @@ std::vector SumRockingEvents(const std::vector & } // namespace -void WriteUnmergedMtzReflections(const std::vector &outcomes, - const UnitCell &unitCell, - const DiffractionExperiment &experiment, - bool sum_partials, - const std::string &filename, - size_t nthreads) { +gemmi::Mtz UnmergedMtz(const std::vector &outcomes, + const UnitCell &unitCell, + const DiffractionExperiment &experiment, + bool sum_partials, + size_t nthreads) { gemmi::Mtz mtz; mtz.spacegroup = &experiment.GetSpaceGroupOrP1(); mtz.set_cell_for_all(unitCell); @@ -885,13 +884,11 @@ void WriteUnmergedMtzReflections(const std::vector &outcomes // stood on that image. Turn the first indexed one back by its own angle to get the orientation of // the sweep, which is the one matrix POINTLESS also writes into every batch. std::optional lattice_at_zero; - std::optional mosaicity_deg; for (const auto &outcome : outcomes) { if (outcome.reflections.empty() || outcome.latt.CalcVolume() <= 1.0f) continue; const float mid_deg = phi_start_deg(outcome.reflections.front().image_number) + wedge_deg / 2.0f; lattice_at_zero = gon ? outcome.latt.Multiply(gon->GetTransformationAngle(mid_deg)) : outcome.latt; - mosaicity_deg = outcome.mosaicity_deg; break; } @@ -930,7 +927,6 @@ void WriteUnmergedMtzReflections(const std::vector &outcomes batch.floats[8 + 3 * i] = u[i] * z_cam; } } - batch.floats[21] = mosaicity_deg.value_or(0.0f); // crydat(0), the reflecting range batch.floats[40] = 1.0f; // scanax = [0, 0, 1]: the rotation axis IS z in the Cambridge frame batch.floats[47] = wedge_deg; batch.floats[61] = 1.0f; // e1 = scanax, the only goniostat axis @@ -975,6 +971,30 @@ void WriteUnmergedMtzReflections(const std::vector &outcomes mtz.data.swap(sorted); mtz.sort_order = {{1, 2, 3, 4, 5}}; } + return mtz; +} + +void SetUnmergedMtzMosaicity(gemmi::Mtz &mtz, const std::vector &outcomes) { + // The reflecting range of the same image the orientation matrix is taken from. + std::optional mosaicity_deg; + for (const auto &outcome : outcomes) { + if (outcome.reflections.empty() || outcome.latt.CalcVolume() <= 1.0f) + continue; + mosaicity_deg = outcome.mosaicity_deg; + break; + } + for (auto &batch : mtz.batches) + batch.floats[21] = mosaicity_deg.value_or(0.0f); // crydat(0), the reflecting range +} + +void WriteUnmergedMtzReflections(const std::vector &outcomes, + const UnitCell &unitCell, + const DiffractionExperiment &experiment, + bool sum_partials, + const std::string &filename, + size_t nthreads) { + gemmi::Mtz mtz = UnmergedMtz(outcomes, unitCell, experiment, sum_partials, nthreads); + SetUnmergedMtzMosaicity(mtz, outcomes); mtz.write_to_file(filename); } diff --git a/image_analysis/WriteReflections.h b/image_analysis/WriteReflections.h index 9aa30e89c..b44722a72 100644 --- a/image_analysis/WriteReflections.h +++ b/image_analysis/WriteReflections.h @@ -11,6 +11,8 @@ #include "../common/DiffractionExperiment.h" #include "IntegrationOutcome.h" +#include + struct MergeStatistics; struct TwinningAnalysisResult; @@ -68,6 +70,16 @@ void WriteUnmergedMtzReflections(const std::vector &outcomes const std::string &filename, size_t nthreads = 1); +// The same file in two steps, for a caller that builds it while the outcomes' per-image scale fields +// are still being written: UnmergedMtz reads everything but mosaicity_deg, which +// SetUnmergedMtzMosaicity then copies into the batch headers. +gemmi::Mtz UnmergedMtz(const std::vector &outcomes, + const UnitCell &unitCell, + const DiffractionExperiment &experiment, + bool sum_partials, + size_t nthreads = 1); +void SetUnmergedMtzMosaicity(gemmi::Mtz &mtz, const std::vector &outcomes); + void WriteReflections(const std::vector &reflections, const UnitCell &unitCell, const DiffractionExperiment &experiment, diff --git a/rugnux/Rugnux.cpp b/rugnux/Rugnux.cpp index 460c308fd..ba9b26063 100644 --- a/rugnux/Rugnux.cpp +++ b/rugnux/Rugnux.cpp @@ -6104,6 +6104,9 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b end_msg.unit_cell = result.consensus_cell; } + // The unmerged MTZ, when it is built beside the P1 cross-check merge (see there). + std::future unmerged_mtz; + // Scaling and merging (full analysis only). if (full && !cancelled_ && (result.indexing_rate.value_or(0.0f) > 0 || result.rotation_lattice_found) && (config_.run_scaling || !config_.reference_data.empty())) { @@ -8742,6 +8745,17 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b const double em_a = result.error_model_a; const double em_b = result.error_model_b; const auto res_fit = result.resolution_fit_A; + // The unmerged MTZ below does not depend on this merge, so it is built meanwhile, from + // the experiment as it stands in the determined group. The merge rewrites each image's + // mosaicity, which the file's batch headers carry, so that is filled in only after it. + if (full && write_files && !geometry_prepass && !superseded && result.consensus_cell + && config_.export_unmerged) + unmerged_mtz = std::async(std::launch::async, + [&outcomes = indexer->GetIntegrationOutcome(), + cell = *result.consensus_cell, x = experiment_, + nthreads = config_.nthreads] { + return UnmergedMtz(outcomes, cell, x, true, nthreads); + }); // Both the merge and the MTZ read the group from the experiment, so it is set for // the whole of it and restored after. experiment_.SpaceGroupNumber(1); @@ -8841,8 +8855,13 @@ ProcessResult Rugnux::RunPipeline(RugnuxObserver *observer, bool write_output, b if (config_.export_unmerged) { if (observer) observer->OnPhase("Writing unmerged reflections"); const std::string path = config_.output_prefix + "_unmerged.mtz"; - WriteUnmergedMtzReflections(indexer->GetIntegrationOutcome(), *result.consensus_cell, - experiment_, true, path, config_.nthreads); + if (unmerged_mtz.valid()) { + gemmi::Mtz mtz = unmerged_mtz.get(); + SetUnmergedMtzMosaicity(mtz, indexer->GetIntegrationOutcome()); + mtz.write_to_file(path); + } else + WriteUnmergedMtzReflections(indexer->GetIntegrationOutcome(), *result.consensus_cell, + experiment_, true, path, config_.nthreads); logger.Info("Unmerged observations written to {}", path); } if (config_.export_unmerged_partials) {