diff --git a/src/aarecommon/math/find_xtal.py b/src/aarecommon/math/find_xtal.py index 1cb0f2d..1640c4d 100644 --- a/src/aarecommon/math/find_xtal.py +++ b/src/aarecommon/math/find_xtal.py @@ -323,9 +323,31 @@ def raster_centre_of_mass(result_array, images) -> CenterOfMassModel | None: return _to_com_model(com, images, "Center of mass") -def raster_highest_score(images) -> CenterOfMassModel | None: - """Target the grid cell with the highest crystal score (compute_crystal_score_array).""" +def raster_highest_score( + images, + min_low_res_spots: float = 10.0, + min_spots_indexed: int = 1, +) -> CenterOfMassModel | None: + """Target the grid cell with the highest crystal score. + + If the whole grid is too weak to hold a crystal — max spots_low_res below + min_low_res_spots OR max spots_indexed below min_spots_indexed — target the + geometric centre of the grid instead of collecting at a noisy cell, so a + 'nothing here' result is centred and deliberate rather than random noise. + """ score_array = compute_crystal_score_array(images) + max_low_res = max((getattr(img, "spots_low_res", 0) or 0 for img in images), default=0) + max_indexed = max((getattr(img, "spots_indexed", 0) or 0 for img in images), default=0) + if max_low_res < min_low_res_spots or max_indexed < min_spots_indexed: + n_nx, n_ny = score_array.shape + logger.info( + f"No crystal in loop (max spots_low_res={max_low_res} < {min_low_res_spots} " + f"or max spots_indexed={max_indexed} < {min_spots_indexed}); " + f"targeting grid centre of {n_nx}x{n_ny} grid" + ) + # get_com_mm maps index i -> (i+0.5)*step, so index (N-1)/2 is the true + # geometric centre of the grid for both odd and even N (and N==1). + return CenterOfMassModel(n_x=(n_nx - 1) / 2.0, n_y=(n_ny - 1) / 2.0) return _to_com_model(_max_cell(score_array), images, "Highest score") @@ -650,3 +672,29 @@ def compare_crystal_methods_scored( ) plt.tight_layout() plt.show() + + +if __name__ == "__main__": + from types import SimpleNamespace as _I + + def _cell(nx, ny, low, idx, bkg=0): + return _I(nx=nx, ny=ny, number=nx * 3 + ny, + spots_low_res=low, spots_indexed=idx, bkg=bkg) + + # Real crystal at (1,1) -> best cell returned. + strong = [_cell(x, y, 0, 0) for x in range(3) for y in range(3)] + strong[4] = _cell(1, 1, 50, 5, bkg=100) + com = raster_highest_score(strong) + assert (round(com.n_x), round(com.n_y)) == (1, 1), com + + # All noise (low-res < 10, nothing indexed) -> centre of 3x3 grid == (1,1). + noise = [_cell(x, y, 3, 0) for x in range(3) for y in range(3)] + com = raster_highest_score(noise) + assert (com.n_x, com.n_y) == (1.0, 1.0), com + + # Strong low-res but nothing indexed anywhere -> OR rule routes to centre. + no_index = [_cell(x, y, 40, 0) for x in range(3) for y in range(3)] + com = raster_highest_score(no_index) + assert (com.n_x, com.n_y) == (1.0, 1.0), com + + print("find_xtal self-check passed")