From c74678b06f415eaeb9db04502b14d0a1c55769c4 Mon Sep 17 00:00:00 2001 From: Filip Leonarski Date: Fri, 23 Oct 2026 16:18:40 +0200 Subject: [PATCH] Polarization plane decided once, in the pre-scan, and followed by the integration too The shadow finder already compares the background's cos 2phi harmonic with the profile a polarized source gives it, and applies the stated degree in the plane the background shows where it shows one clearly (correlation > 0.5, at least 0.3 of a fully polarized source's modulation). The integration's Lorentz-polarization factor still assumed the plane along lab x, so on a run whose background is suppressed along the image's y axis the two disagreed, and only a negative --polarization made them agree. ShadowFinder now reports what it applied (GetAppliedPolarization), the pre-scan writes it into the experiment (and keeps it in the beam-stop memo with the mask), and everything built after the pre-scan - the Bragg integration, the azimuthal mapping - reads the same plane. The degree is never fitted; only the plane can turn, and only where the shadow finder already turned it. Where the plane is not turned the experiment is not touched, so every other run is bit-identical: six rotation sets give the same p.hkl md5 as rc178, and the stills sets that do not turn give the same R_free, R_split, indexed fraction and d_min. On the one stills run whose plane turns, R_free and the resolution cut are the ones a stated -0.998 gave before, now from a bare invocation. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_015CeEGgpKdfkCTJb3YNMJcu (cherry picked from commit f50bfe15f7dccf588bd8bdd6114a32e3ac9b36ed) --- docs/CPU_DATA_ANALYSIS_IMAGE.md | 7 ++++++- docs/RUGNUX_ADVANCED.md | 2 +- image_analysis/beam_stop/SHADOW_FINDER.md | 4 +++- image_analysis/beam_stop/ShadowFinder.cpp | 11 ++++++++++- image_analysis/beam_stop/ShadowFinder.h | 7 +++++++ rugnux/Rugnux.h | 1 + rugnux/RugnuxPreScan.cpp | 15 ++++++++++++++- rugnux/rugnux_cli.cpp | 2 +- tests/ShadowFinderTest.cpp | 12 +++++++++++- 9 files changed, 54 insertions(+), 7 deletions(-) diff --git a/docs/CPU_DATA_ANALYSIS_IMAGE.md b/docs/CPU_DATA_ANALYSIS_IMAGE.md index 298bce3bc..45423a635 100644 --- a/docs/CPU_DATA_ANALYSIS_IMAGE.md +++ b/docs/CPU_DATA_ANALYSIS_IMAGE.md @@ -495,7 +495,12 @@ determine it, so it comes from the file where the file states it: a miniCBF `# P a d*TREK `SOURCE_POLARZ` both give the fraction $f$, read as $P = 2f - 1$ (so a laboratory tube's 0.5 is unpolarized), and `SOURCE_POLARZ` also gives the plane's normal, which sets the sign. Where the file says nothing (ADSC SMV, marCCD, NXmx) the default $P = 0.99$, an undulator value, applies, with the -plane horizontal; a vertical beamline must say so with `--polarization`. Bending-magnet and wiggler +plane horizontal. The **plane**, unlike the degree, is measured: the beam-stop pre-scan compares the +background's $\cos 2\phi$ harmonic with the profile a polarized source gives it, and where the background +clearly follows the other plane (correlation above 0.5 and at least 0.3 of a fully polarized source's +modulation) the stated degree is applied in that plane - to the shadow test and, from then on, to the +integration and every other correction of the run. Elsewhere the plane stays as stated, so a run whose +beam stop is not detected (`--detect-beam-stop=off`) keeps the stated one. Bending-magnet and wiggler beamlines are typically lower (0.8–0.95), and a known value for one should be passed with `--polarization`. diff --git a/docs/RUGNUX_ADVANCED.md b/docs/RUGNUX_ADVANCED.md index cfa49d21f..308f6a945 100644 --- a/docs/RUGNUX_ADVANCED.md +++ b/docs/RUGNUX_ADVANCED.md @@ -504,5 +504,5 @@ Geometry overrides (defaults are taken from the input file; override them to rep | `--rot3 ` | PONI detector rotation 3, about the beam (rad) — the one `rugnux --mode calibration` writes into its `.poni`, so a refined third rotation can be given straight back | | `--detector-mirror-y` | The stored image is mirrored in Y relative to the frame the PONI angles are stated in | | `--detector-quarter-turns <0-3>` | The stored image is turned by this many multiples of 90° about the beam relative to that same frame | -| `--polarization ` | Degree of polarization *p* (XDS `FRACTION_OF_POLARIZATION` = (1+*p*)/2), signed: positive for a horizontally polarized source, negative for a vertically polarized one. Default: what a miniCBF `# Polarization` line or a d*TREK `SOURCE_POLARZ` states (a fraction *f*, read as *p* = 2*f* − 1, so a laboratory tube's 0.5 is unpolarized), and otherwise 0.99, an undulator value; bending-magnet and wiggler sources are typically lower (0.8–0.95), and neither ADSC SMV, marCCD nor NXmx files state it, so pass the beamline's known value here | +| `--polarization ` | Degree of polarization *p* (XDS `FRACTION_OF_POLARIZATION` = (1+*p*)/2), signed: positive for a horizontally polarized source, negative for a vertically polarized one; where the beam-stop pre-scan finds the background clearly suppressed in the other plane, the degree is applied in that plane instead. Default: what a miniCBF `# Polarization` line or a d*TREK `SOURCE_POLARZ` states (a fraction *f*, read as *p* = 2*f* − 1, so a laboratory tube's 0.5 is unpolarized), and otherwise 0.99, an undulator value; bending-magnet and wiggler sources are typically lower (0.8–0.95), and neither ADSC SMV, marCCD nor NXmx files state it, so pass the beamline's known value here | | `--rotation-scale ` | Goniometer rotation scale: the stage turned `k` times the angle stored in the file (the commanded one). Applied to both passes, and overrides the scale Rugnux fits for itself | diff --git a/image_analysis/beam_stop/SHADOW_FINDER.md b/image_analysis/beam_stop/SHADOW_FINDER.md index 998947ea8..fd451c6e2 100644 --- a/image_analysis/beam_stop/SHADOW_FINDER.md +++ b/image_analysis/beam_stop/SHADOW_FINDER.md @@ -213,7 +213,9 @@ beam-stop mask is an **offline** (Rugnux + viewer) feature and never goes to the above 0.5 and at least 0.3 of a fully polarized source's modulation - and as stated otherwise. Over the rotation and stills corpora only the run whose plane lies along the image's y axis is turned; on every other run the correction, and so the mask, is what it was. An experiment that states no - polarization is left uncorrected, as before. + polarization is left uncorrected, as before. What was applied is read back with + `GetAppliedPolarization()`, and `rugnux`'s pre-scan writes it into the experiment, so the + integration's correction follows the same plane: the plane is decided once, here. - **Unsigned saturation**: for unsigned images a detector-saturating reflection collides with the no-data sentinel and is skipped; harmless in practice (such reflections do not sit under the stop), and absent entirely on the signed/preprocessed path. diff --git a/image_analysis/beam_stop/ShadowFinder.cpp b/image_analysis/beam_stop/ShadowFinder.cpp index a042a40c6..171a569ba 100644 --- a/image_analysis/beam_stop/ShadowFinder.cpp +++ b/image_analysis/beam_stop/ShadowFinder.cpp @@ -707,6 +707,11 @@ void ShadowFinder::ReleaseProjection() { #endif } +std::optional ShadowFinder::GetAppliedPolarization() const { + std::unique_lock ul(m); + return applied_polarization; +} + std::vector ShadowFinder::GetMeanProjection() const { std::unique_lock ul(m); #ifdef JFJOCH_USE_CUDA @@ -732,6 +737,7 @@ std::vector ShadowFinder::GetMask(size_t nthreads) const { std::unique_lock ul(m); if (nthreads == 0) nthreads = AvailableCpus(); + applied_polarization.reset(); #ifdef JFJOCH_USE_CUDA // Where every frame went to the device the mask is made there too, from the projection as it lies - // when it lies as the detector does (no panel table). @@ -752,11 +758,13 @@ std::vector ShadowFinder::GetMask(size_t nthreads) const { setup.distance_mm = geometry.GetDetectorDistance_mm(); ScopedGpuPin pin(gpu->GetDevice()); setup.has_polarization = polarization.has_value(); - if (polarization) + if (polarization) { setup.polarization = PlaneApplied(*polarization, measured_plane(gpu->MeanProjection(pixel_mask), width, height, beam_x, beam_y, setup.det_matrix, setup.distance_mm / setup.pixel_size_mm, nthreads)); + applied_polarization = setup.polarization; + } return gpu->Mask(setup, pixel_mask); } #endif @@ -825,6 +833,7 @@ std::vector ShadowFinder::MaskOn(const Projection &p, const std::vecto }); polarization_applied = PlaneApplied(*polarization, measured_plane(projection, W, H, centre_x, centre_y, det_matrix.data(), distance_px, nthreads)); + applied_polarization = polarization_applied; } Plane mean(n_pixels); diff --git a/image_analysis/beam_stop/ShadowFinder.h b/image_analysis/beam_stop/ShadowFinder.h index 1c844f5cd..d87e59bf5 100644 --- a/image_analysis/beam_stop/ShadowFinder.h +++ b/image_analysis/beam_stop/ShadowFinder.h @@ -61,6 +61,8 @@ class ShadowFinder { // it lies. Only the centre is replaced, by the one the caller measured. const DiffractionGeometry geometry; const std::optional polarization; // unset leaves the background as it was measured + // What the last GetMask divided by: the stated degree, in the plane the background shows. + mutable std::optional applied_polarization; std::vector pixel_mask; // pixels already masked carry no background to test // The panels the image is made of. With any, the mask is made on the detector plane (GetMask). @@ -156,6 +158,11 @@ public: // dominate this, and they are all exactly parallel. [[nodiscard]] std::vector GetMask(size_t nthreads = 0) const; + // The polarization the last GetMask divided the background by: the stated degree, in the plane the + // background shows where it shows one clearly, and as stated otherwise. The run's other corrections + // take it from here, so that the plane is decided once. Unset with no stated degree, or before GetMask. + [[nodiscard]] std::optional GetAppliedPolarization() const; + // Mean counts per pixel over the frames added, NAN where nothing was counted. This is the // projection GetMask() tests, so anything else that wants the background before indexing // gets it without reading the frames a second time. diff --git a/rugnux/Rugnux.h b/rugnux/Rugnux.h index 9b539dfc2..962558a1a 100644 --- a/rugnux/Rugnux.h +++ b/rugnux/Rugnux.h @@ -1383,6 +1383,7 @@ class Rugnux { std::string centre_line; // what the detection said about the centre size_t frames_read = 0; size_t shadow_frames = 0; + std::optional polarization; // ShadowFinder::GetAppliedPolarization }; std::mutex m; std::vector entries; diff --git a/rugnux/RugnuxPreScan.cpp b/rugnux/RugnuxPreScan.cpp index cc2154361..e84061b90 100644 --- a/rugnux/RugnuxPreScan.cpp +++ b/rugnux/RugnuxPreScan.cpp @@ -1011,9 +1011,11 @@ void Rugnux::PreScan(int start_image, int images_to_process, int frame_count, Ru // consumes is the one fitted below instead, with the shadow out of the way. std::vector shadow; size_t shadow_frames = 0; + std::optional polarization; if (shadow_memo) { projection = shadow_memo->projection; ring_centre = shadow_memo->ring_centre; + polarization = shadow_memo->polarization; shadow.assign(experiment_.GetXPixelsNumConv() * experiment_.GetYPixelsNumConv(), 0); for (const auto &[i, v] : shadow_memo->shadow) shadow[i] = v; @@ -1041,10 +1043,11 @@ void Rugnux::PreScan(int start_image, int images_to_process, int frame_count, Ru logger.Info("{}", centre_line); shadow = finder->GetMask(config_.nthreads); shadow_frames = finder->GetFrameCount(); + polarization = finder->GetAppliedPolarization(); finder.reset(); BeamStopMemo::Entry entry{shadow_key, {}, projection, ring_centre, centre_line, frames_read.load(), - shadow_frames}; + shadow_frames, polarization}; for (size_t i = 0; i < shadow.size(); i++) if (shadow[i] != 0) entry.shadow.emplace_back(static_cast(i), shadow[i]); @@ -1055,6 +1058,16 @@ void Rugnux::PreScan(int start_image, int images_to_process, int frame_count, Ru entries.erase(entries.begin()); entries.push_back(std::move(entry)); } + // The polarization plane, decided once: the shadow finder has just compared the background with + // the stated plane and turned it where the background clearly lies the other way (ShadowFinder), + // and the integration's correction follows the same plane from here on. + if (polarization && polarization != experiment_.GetPolarizationFactor()) { + logger.Info("Polarization: the background is suppressed in the plane turned by 90 degrees from the " + "stated one, so the correction is applied in that plane (coefficient {:.3f})", + *polarization); + Note("polarization plane taken from the background, turned by 90 degrees"); + experiment_.PolarizationFactor(polarization); + } const auto shadowed = std::count_if(shadow.begin(), shadow.end(), [](uint32_t v) { return v != 0; }); pixel_mask_.LoadBeamStopMask(experiment_, shadow); opaque_shadow = shadow; diff --git a/rugnux/rugnux_cli.cpp b/rugnux/rugnux_cli.cpp index 9a9f45e01..3ef431b71 100644 --- a/rugnux/rugnux_cli.cpp +++ b/rugnux/rugnux_cli.cpp @@ -284,7 +284,7 @@ void print_usage() { std::cout << " --rot3 PONI rotation 3, about the beam (rad)" << std::endl; std::cout << " --detector-mirror-y Stored image is mirrored in Y vs the detector" << std::endl; std::cout << " --detector-quarter-turns <0-3> Stored image is turned by this many 90 deg about the beam" << std::endl; - std::cout << " --polarization Degree of polarization p (XDS FRACTION_OF_POLARIZATION = (1+p)/2), signed: positive for a horizontally polarized source, negative for a vertically polarized one (the sign rotates the polarization plane by 90 degrees). Default: what a miniCBF \"# Polarization\" line or a d*TREK SOURCE_POLARZ states (a fraction f, read as p = 2f - 1, so 0.5 from a laboratory tube is unpolarized), and otherwise 0.99, an undulator value; bending-magnet and wiggler sources are typically lower (0.8-0.95), and neither ADSC SMV, marCCD nor NXmx files state it, so pass a known beamline value here" << std::endl; + std::cout << " --polarization Degree of polarization p (XDS FRACTION_OF_POLARIZATION = (1+p)/2), signed: positive for a horizontally polarized source, negative for a vertically polarized one (the sign rotates the polarization plane by 90 degrees); where the beam-stop pre-scan finds the background clearly suppressed in the other plane, the degree is applied in that plane instead. Default: what a miniCBF \"# Polarization\" line or a d*TREK SOURCE_POLARZ states (a fraction f, read as p = 2f - 1, so 0.5 from a laboratory tube is unpolarized), and otherwise 0.99, an undulator value; bending-magnet and wiggler sources are typically lower (0.8-0.95), and neither ADSC SMV, marCCD nor NXmx files state it, so pass a known beamline value here" << std::endl; } enum { diff --git a/tests/ShadowFinderTest.cpp b/tests/ShadowFinderTest.cpp index 8966e4e43..4a4b9384a 100644 --- a/tests/ShadowFinderTest.cpp +++ b/tests/ShadowFinderTest.cpp @@ -412,6 +412,7 @@ TEST_CASE("ShadowFinder_APolarizedBackgroundIsNotAShadow", "[ShadowFinder]") { frame[I(x, y)] = static_cast(std::lround(LEVEL * factor / at_peak)); } + std::optional applied; auto mask_for = [&](const std::optional &polarization) { DiffractionExperiment x = TestExperiment(); x.DetectorDistance_mm(DISTANCE_MM).PolarizationFactor(polarization); @@ -424,10 +425,12 @@ TEST_CASE("ShadowFinder_APolarizedBackgroundIsNotAShadow", "[ShadowFinder]") { finder.AddImage(msg, buffer); } const auto mask = finder.GetMask(); + applied = finder.GetAppliedPolarization(); return std::count(mask.begin(), mask.end(), 1u); }; CHECK(mask_for(POLARIZATION) == 0); + CHECK(applied == POLARIZATION); // the plane the background shows is the stated one: kept as stated CHECK(mask_for(std::nullopt) > static_cast(W) * H / 20); } @@ -455,6 +458,7 @@ TEST_CASE("ShadowFinder_ThePolarizationPlaneIsTheOneTheBackgroundShows", "[Shado frame[I(xi, y)] = stop[I(xi, y)] == 0 ? 0 : static_cast(std::lround(LEVEL * factor)); } + std::optional applied; // what the run's other corrections take (GetAppliedPolarization) auto mask_for = [&](const std::optional &polarization) { DiffractionExperiment x = TestExperiment(); x.DetectorDistance_mm(DISTANCE_MM).PolarizationFactor(polarization); @@ -466,10 +470,13 @@ TEST_CASE("ShadowFinder_ThePolarizationPlaneIsTheOneTheBackgroundShows", "[Shado msg.image = CompressedImage(frame, W, H); finder.AddImage(msg, buffer); } - return finder.GetMask(); + const auto out = finder.GetMask(); + applied = finder.GetAppliedPolarization(); + return out; }; const auto mask = mask_for(0.99f); + CHECK(applied == -0.99f); // the stated degree, in the plane shown CHECK(mask[I(C, C)] == ShadowFinder::SHADOW); // the stop CHECK(mask[I(W - 3, C)] == ShadowFinder::SHADOW); // and its arms, along x CHECK(mask[I(C, 3)] == ShadowFinder::SHADOW); // and along y, inside the suppressed lobe @@ -477,6 +484,9 @@ TEST_CASE("ShadowFinder_ThePolarizationPlaneIsTheOneTheBackgroundShows", "[Shado CHECK(mask[I(C - 20, H - 4)] == 0); CHECK(mask_for(0.998f) == mask); CHECK(mask_for(-0.998f) == mask); + CHECK(applied == -0.998f); + mask_for(std::nullopt); + CHECK(!applied); // no stated degree, nothing to apply } // The plane is decided from the sector medians of each band: the cos 2phi harmonic has to follow the