diff --git a/src/aare/common/find_xtal.py b/src/aare/common/find_xtal.py index c149434c..0e3acf0e 100644 --- a/src/aare/common/find_xtal.py +++ b/src/aare/common/find_xtal.py @@ -280,32 +280,35 @@ def get_com_image_number(com, images): return image.number return None -def raster_centre_of_mass(result_array, images) -> CenterOfMassModel | None: - # grid_mm_x and grid_mm_y are relative to the top left corner of raster grid - # com = ndimage.center_of_mass(result_array) - # if com and not np.isnan(com[0]) and not np.isnan(com[1]): - # logger.info(f"Center of mass: {com}") - # max_image = get_com_image_number(com, images) - # return CenterOfMassModel(n_x=com[0], n_y=com[1], max_image=max_image) - # elif com and np.isnan(com[0]) and np.isnan(com[1]): - # logger.warning(f"Center of mass is nan: {com[0]}, {com[1]}") - # return None - # else: - # logger.warning("No valid center of mass found") - # return None - score_array = compute_crystal_score_array(images) - com = _max_cell(score_array) - if com and not np.isnan(com[0]) and not np.isnan(com[1]): - logger.info(f"Highest score: {com}") - max_image = get_com_image_number(com, images) - return CenterOfMassModel(n_x=com[0], n_y=com[1], max_image=max_image) - elif com and np.isnan(com[0]) and np.isnan(com[1]): - logger.warning(f"Highest score is nan: {com[0]}, {com[1]}") +def _to_com_model(coords, images, label: str) -> CenterOfMassModel | None: + """Validate (n_x, n_y) grid coords and wrap them in a CenterOfMassModel. + + Shared by raster_centre_of_mass and raster_highest_score, which differ only + in how they pick the target cell. + """ + if coords and not np.isnan(coords[0]) and not np.isnan(coords[1]): + logger.info(f"{label}: {coords}") + max_image = get_com_image_number(coords, images) + return CenterOfMassModel(n_x=coords[0], n_y=coords[1], max_image=max_image) + elif coords and np.isnan(coords[0]) and np.isnan(coords[1]): + logger.warning(f"{label} is nan: {coords[0]}, {coords[1]}") return None else: - logger.warning("No valid score found") + logger.warning(f"No valid {label} found") return None + +def raster_centre_of_mass(result_array, images) -> CenterOfMassModel | None: + # grid_mm_x and grid_mm_y are relative to the top left corner of raster grid + com = ndimage.center_of_mass(result_array) + 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).""" + score_array = compute_crystal_score_array(images) + return _to_com_model(_max_cell(score_array), images, "Highest score") + def has_sufficient_low_res_spots( result_array: np.ndarray, min_spots_low_res: float, diff --git a/src/aare/common/raster_grid.py b/src/aare/common/raster_grid.py index 082141c5..67c6c559 100644 --- a/src/aare/common/raster_grid.py +++ b/src/aare/common/raster_grid.py @@ -82,6 +82,7 @@ class RasterPayloadModel(BaseModel): sample_id: int attach_image: bool = True centre_of_mass: CenterOfMassModel | None = None + raster_score: List[float | None] | None = None center_pxl: Coordinate | None start_pxl: Coordinate cell_size_pxl: Annotated[Coordinate, AfterValidator(positive_coords)] diff --git a/src/aare/daq/aaredb.py b/src/aare/daq/aaredb.py index 035aa750..472616e4 100644 --- a/src/aare/daq/aaredb.py +++ b/src/aare/daq/aaredb.py @@ -25,6 +25,7 @@ from aareDB import ( from pydantic import StrictInt from aare.common.coordinate import Coordinate +from aare.common.find_xtal import compute_crystal_score_array from aare.common.logger_config import setup_logger from aare.common.logger_events import log_timing from aare.common.models import ( @@ -358,12 +359,21 @@ class AareWrapper: else: center_pxl = None + try: + score_arr = compute_crystal_score_array(raster_result.images) + score = [float(score_arr[img.nx, img.ny]) if img.nx is not None and img.ny is not None else None + for img in raster_result.images] + except Exception as e: + logger.warning(f"raster score computation failed, sending null score: {e}") + score = None + payload = RasterPayloadModel( request = raster_request, result = raster_result, sample_id = sample.db_id, attach_image = True, centre_of_mass = com, + raster_score = score, # TODO: Check whether AareDB needs to change for taking this input start_pxl = start_pxl, center_pxl = center_pxl, cell_size_pxl = cell_size_pxl, diff --git a/src/aare/daq/operations/raster/service.py b/src/aare/daq/operations/raster/service.py index 8b5804b3..0c361ef1 100644 --- a/src/aare/daq/operations/raster/service.py +++ b/src/aare/daq/operations/raster/service.py @@ -8,7 +8,7 @@ from jfjoch_client.exceptions import NotFoundException from aare.common.coordinate import AerotechCoordinate, Coordinate, SmargonCoordinate from aare.common.exception_handler import AutoRasterSampleSkipped, RasterScanException -from aare.common.find_xtal import create_quality_filtered_array, raster_centre_of_mass +from aare.common.find_xtal import raster_highest_score from aare.common.logger_events import ( geom_log_context, log_ml_bundle_meta, @@ -351,16 +351,7 @@ class RasterService: if self.ctx.deps.cfg.simulated_detector: self.logger.info("Simulated detector mode enabled; using fake raster result.") scan_result = generate_no_beam_scan_result(request) - result_array = create_quality_filtered_array( - scan_result.images, - "spots_low_res", - min_spots=None, - min_efficiency=1.0, - min_background=None, - min_low_res_spots=None, - ) - - com = raster_centre_of_mass(result_array, scan_result.images) + com = raster_highest_score(scan_result.images) target_coor = com.get_com_mm(request) target_coor_offset = self.ctx.sample_geometry.smargon_nudge(target_coor) target_smargon = SmargonCoordinate( @@ -381,15 +372,7 @@ class RasterService: ) raise RasterScanException("JFJoch returned no ScanResult for raster") - result_array = create_quality_filtered_array( - scan_result.images, - "spots_low_res", - min_spots=None, - min_efficiency=1.0, - min_background=None, - min_low_res_spots=None, - ) - com = raster_centre_of_mass(result_array, scan_result.images) + com = raster_highest_score(scan_result.images) if com is None: self.logger.info("Calcualted COM is None using centre image") if request.n_x == 1: diff --git a/src/aare/gui/panels/raster_data_collection.py b/src/aare/gui/panels/raster_data_collection.py index fc6f65a9..ce214bdd 100644 --- a/src/aare/gui/panels/raster_data_collection.py +++ b/src/aare/gui/panels/raster_data_collection.py @@ -86,6 +86,7 @@ class RasterDataCollectionPanel(ScanSettingsPanel): self.metric_combo = QComboBox() self._layout.addWidget(QLabel("Metric", parent=self), 7, 0) + self.metric_combo.addItem("Raster score", RasterGridMetric.RASTER_SCORE) self.metric_combo.addItem("Spot count (low res.)", RasterGridMetric.SPOTS_LOW_RES) self.metric_combo.addItem("Spot count", RasterGridMetric.SPOTS) self.metric_combo.addItem("Spot count (indexed)", RasterGridMetric.SPOTS_INDEXED) @@ -93,11 +94,11 @@ class RasterDataCollectionPanel(ScanSettingsPanel): self.metric_combo.addItem("Spot ratio (ice/low res.)", RasterGridMetric.SPOTS_ICE_LOW_RES) self.metric_combo.addItem("Background estimate", RasterGridMetric.BKG) self.metric_combo.addItem("Indexing result", RasterGridMetric.INDEXING) - self.metric_combo.addItem("Profile Radius (Score)", RasterGridMetric.PR) # JD: repurposing PR for testing scoring + self.metric_combo.addItem("Profile Radius", RasterGridMetric.PR) self.metric_combo.addItem("B-factor", RasterGridMetric.BFACTOR) self.metric_combo.addItem("Resolution (ML est.)", RasterGridMetric.RES) - self.metric_combo.setCurrentIndex(self.metric_combo.findData(RasterGridMetric.PR)) # JD: repurposing PR for testing scoring + self.metric_combo.setCurrentIndex(self.metric_combo.findData(RasterGridMetric.RASTER_SCORE)) self.metric_combo.currentIndexChanged.connect(self.metric_changed) self._layout.addWidget(self.metric_combo, 7, 1, 1, 3) diff --git a/src/aare/gui/scan_logic/raster_grid_manager.py b/src/aare/gui/scan_logic/raster_grid_manager.py index ca1503da..122e6abf 100644 --- a/src/aare/gui/scan_logic/raster_grid_manager.py +++ b/src/aare/gui/scan_logic/raster_grid_manager.py @@ -30,6 +30,7 @@ class RasterGridMetric(Enum): SPOTS_ICE = 7 SPOTS_INDEXED = 8 SPOTS_ICE_LOW_RES = 9 + RASTER_SCORE = 10 # Viridis colormap colors (from matplotlib), dark purple -> yellow _VIRIDIS_COLORS = [ @@ -122,8 +123,7 @@ class RasterGridManager(QObject): self.__beam_size_mm = geom.beam_size_mm self.__geom = geom - # self.__metric = RasterGridMetric.SPOTS_LOW_RES - self.__metric = RasterGridMetric.PR # JD: repurposing PR for testing scoring + self.__metric = RasterGridMetric.RASTER_SCORE self.__loaded_image_prefix = None self.__loaded_image_index = None @@ -375,11 +375,11 @@ class RasterGridManager(QObject): if grid.result.images[cell].b is not None and grid.result.images[cell].b >= 0: txt += f"B-factor {grid.result.images[cell].b:.2f}
" - # if grid.result.images[cell].pr is not None and grid.result.images[cell].pr >= 0: # JD: repurposing PR for testing scoring - # txt += f"Profile Radius {grid.result.images[cell].pr:.2f}" + if grid.result.images[cell].pr is not None and grid.result.images[cell].pr >= 0: + txt += f"Profile Radius {grid.result.images[cell].pr:.2f}
" if grid.result.images[cell].nx is not None and grid.result.images[cell].ny is not None: score = compute_crystal_score_array(grid.result.images) - txt += f"Profile Radius (Score) {score[grid.result.images[cell].nx, grid.result.images[cell].ny]:.2f}" + txt += f"Raster score {score[grid.result.images[cell].nx, grid.result.images[cell].ny]:.2f}" return txt return None @@ -466,8 +466,9 @@ class RasterGridManager(QObject): v = [obj.bkg for obj in i.result.images] case RasterGridMetric.INDEXING: v = [obj.index for obj in i.result.images] - case RasterGridMetric.PR: # JD: repurposing PR for testing scoring - # v = [obj.index / max(obj.spots_low_res, 1) for obj in i.result.images] + case RasterGridMetric.PR: + v = [obj.index / max(obj.spots_low_res, 1) for obj in i.result.images] + case RasterGridMetric.RASTER_SCORE: score = compute_crystal_score_array(i.result.images) v = [score[obj.nx, obj.ny] if obj.nx is not None and obj.ny is not None else None for obj in i.result.images] case RasterGridMetric.BFACTOR: diff --git a/tests/unit/common/test_find_xtal.py b/tests/unit/common/test_find_xtal.py index 13b74e96..d7cd264c 100644 --- a/tests/unit/common/test_find_xtal.py +++ b/tests/unit/common/test_find_xtal.py @@ -4,7 +4,8 @@ from unittest.mock import MagicMock from aare.common.find_xtal import ( identify_crystal_raster, rebuild_array_from_scan_results, create_quality_filtered_array, get_xtal_size, get_best_b_factor, - get_best_res, raster_centre_of_mass, has_sufficient_low_res_spots + get_best_res, raster_centre_of_mass, has_sufficient_low_res_spots, + compute_crystal_score_array, raster_highest_score ) from aare.common.raster_grid import RasterGridRequest, CenterOfMassModel from aare.common.models import CrystalSize, Coordinate @@ -95,3 +96,51 @@ def test_has_sufficient_low_res_spots(): assert has_sufficient_low_res_spots(arr, 3.0) is True assert has_sufficient_low_res_spots(arr, 5.0) is False assert has_sufficient_low_res_spots(None, 1.0) is False + +@pytest.fixture +def mock_score_results(): + # 2x2 grid; cell (1, 1) is the strongest crystal across bkg, spots_low_res + # and spots_indexed (the three inputs to compute_crystal_score_array). + data = [ + # nx, ny, bkg, spots_low_res, spots_indexed + (0, 0, 1.0, 1.0, 0.0), + (1, 0, 2.0, 2.0, 0.0), + (0, 1, 2.0, 2.0, 0.0), + (1, 1, 10.0, 20.0, 5.0), + ] + results = [] + for i, (nx, ny, bkg, low, idx) in enumerate(data): + res = MagicMock() + res.nx, res.ny, res.number = nx, ny, i + res.bkg, res.spots_low_res, res.spots_indexed = bkg, low, idx + results.append(res) + return results + +def test_compute_crystal_score_array(mock_score_results): + score = compute_crystal_score_array(mock_score_results) + assert score.shape == (2, 2) + # weights sum to 1.0 and each input is min-max normalised to [0, 100] + assert score.min() >= 0.0 and score.max() <= 100.0 + # (1, 1) is max in all three inputs -> 100; weak corner (0, 0) is min in all -> 0 + assert score[1, 1] == pytest.approx(100.0) + assert score[0, 0] == pytest.approx(0.0) + assert np.unravel_index(np.argmax(score), score.shape) == (1, 1) + +def test_compute_crystal_score_array_weights(): + # Cell A is the sole max in spots_low_res (weight 0.55); cell B is the sole + # max in spots_indexed (weight 0.20). A must outscore B. + def _cell(nx, n, low, idx): + r = MagicMock() + r.nx, r.ny, r.number = nx, 0, n + r.bkg, r.spots_low_res, r.spots_indexed = 0.0, low, idx + return r + score = compute_crystal_score_array([_cell(0, 0, 10.0, 0.0), _cell(1, 1, 0.0, 10.0)]) + assert score[0, 0] > score[1, 0] + assert score[0, 0] == pytest.approx(55.0) + assert score[1, 0] == pytest.approx(20.0) + +def test_raster_highest_score(mock_score_results): + com = raster_highest_score(mock_score_results) + assert com.n_x == 1.0 + assert com.n_y == 1.0 + assert com.max_image == 3