BeamCenterFromBackground: skip pixels clearly outside the band before sqrt and atan2

60-70% of a large detector lies outside the resolution band the walk bins, and every one paid a
square root and two arctangents on every iteration. A test on tan(2theta) = rho / lz with a 0.1%
margin, far above float rounding, drops them first; every pixel the exact test keeps still reaches it.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01D1G8gJVAy6gp1K5Dz3NE5C

image_analysis/geom_refinement/BeamCenterFromBackground.cpp

@@ -127,6 +127,12 @@ FindBeamCenterFromBackground(const DiffractionExperiment &experiment, const Pixe
     std::vector<double> block_jy(static_cast<size_t>(BLOCKS) * n_cells);
     std::vector<int32_t> block_count(static_cast<size_t>(BLOCKS) * n_cells);
 
+    // Most pixels lie outside the band. Those clearly outside it in tan(2theta) = rho / lz - by a
+    // margin far above float rounding - skip the square root and both atan2 below; every pixel the exact
+    // test would keep still reaches it.
+    const double tan_lo = std::tan(static_cast<double>(tt_lo)) * (1.0 - 1e-3);
+    const double tan_hi = tt_hi < PI / 2 ? std::tan(static_cast<double>(tt_hi)) * (1.0 + 1e-3) : INFINITY;
+
     float step_x = 0.0f, step_y = 0.0f, sigma_x = 0.0f, sigma_y = 0.0f;
     float previous_x = 0.0f, previous_y = 0.0f;
     int reversals = 0;
@@ -153,7 +159,13 @@ FindBeamCenterFromBackground(const DiffractionExperiment &experiment, const Pixe
                     const float lx = rot[0] * u + rot[1] * v + rot[2] * distance;
                     const float ly = rot[3] * u + rot[4] * v + rot[5] * distance;
                     const float lz = rot[6] * u + rot[7] * v + rot[8] * distance;
-                    const float rho = std::sqrt(lx * lx + ly * ly);
+                    const float rho_sq = lx * lx + ly * ly;
+                    if (lz > 0.0f) {
+                        const double lz_sq = static_cast<double>(lz) * lz;
+                        if (rho_sq < tan_lo * tan_lo * lz_sq || rho_sq > tan_hi * tan_hi * lz_sq)
+                            continue;
+                    }
+                    const float rho = std::sqrt(rho_sq);
                     const float two_theta = std::atan2(rho, lz);
                     if (two_theta < tt_lo || two_theta >= tt_hi || rho == 0.0f)
                         continue;