Fix openslide offset

src/classpose/entrypoints/predict_wsi.py

@@ -146,6 +146,8 @@ def __init__(
         self.ts = manager.Value("f", 0)
         self.mpp_x = manager.Value("f", 0)
         self.mpp_y = manager.Value("f", 0)
+        self.bounds_x = manager.Value("f", 0)
+        self.bounds_y = manager.Value("f", 0)
         self.tissue_cnts = manager.list()
         self.roi_cnts = manager.list()
 
@@ -169,6 +171,17 @@ def _init_slide(self):
         self.mpp = get_slide_resolution(self.slide)
         self.mpp_x.value = self.mpp[0]
         self.mpp_y.value = self.mpp[1]
+        bounds_x = self.slide.properties.get("openslide.bounds-x")
+        bounds_y = self.slide.properties.get("openslide.bounds-y")
+        self.bounds_x.value = float(bounds_x) if bounds_x is not None else 0.0
+        self.bounds_y.value = float(bounds_y) if bounds_y is not None else 0.0
+        logger.info(f"Slide bounds offset: x={self.bounds_x.value}, y={self.bounds_y.value}")
+
+        if self.roi_tree is not None and (
+            self.bounds_x.value != 0 or self.bounds_y.value != 0
+        ):
+            self._align_roi_tree_to_slide_bounds()
+
         target_downsample = min(
             self.train_mpp / self.mpp[0], self.train_mpp / self.mpp[1]
         )
@@ -196,6 +209,27 @@ def _init_slide(self):
         logger.info(f"Slide downsample: {self.ts.value}")
         logger.info(f"Level: {self.level}")
 
+    def _align_roi_tree_to_slide_bounds(self):
+        """
+        Align ROI polygons to slide coordinates when OpenSlide bounds offsets are present.
+
+        QuPath ROI annotations are exported in coordinates relative to the displayed image 
+        origin (i.e. bounds-offset corrected). Internal OpenSlide tile coordinates are in 
+        level-0 slide space. So if needed, shift ROI polygons by +bounds offsets so 
+        ROI-based tile selection/filtering uses the same frame.
+        """
+        bounds_x_val = float(self.bounds_x.value)
+        bounds_y_val = float(self.bounds_y.value)
+        geometries = list(self.roi_tree.geometries)
+
+        logger.info(f"Applying bounds offset to ROI polygons: x={bounds_x_val}, y={bounds_y_val}")
+
+        shifted = [
+            shapely.affinity.translate(geom, xoff=bounds_x_val, yoff=bounds_y_val)
+            for geom in geometries
+        ]
+        self.roi_tree = shapely.STRtree(shifted)
+
     def _get_tissue_contours(self):
         if self.tissue_detection_model_path is not None:
             logger.info("Detecting tissue contours using GrandQC")
@@ -460,7 +494,7 @@ def __call__(
                     continue
                 center = np.round(polygon.centroid.coords[0], 2).tolist()
                 curr_coords = curr_coords.tolist()
-                curr_coords.append(curr_coords[0])
+                curr_coords.append(curr_coords[0].copy())
                 cl = class_masks[cell_mask][0]
                 curr_cell = {
                     "id": str(uuid.uuid4()),
@@ -611,6 +645,44 @@ def to_geojson_polygon(curr_cell: dict) -> dict:
     return curr_cell
 
 
+def apply_bounds_offset_to_feature(feature: dict, bounds_x: float, bounds_y: float) -> dict:
+    """
+    Apply bounds offset to a GeoJSON Polygon feature's geometry and centroid measurements.
+
+    Args:
+        feature: GeoJSON Feature dict with Polygon geometry
+        bounds_x: X offset from openslide.bounds-x property
+        bounds_y: Y offset from openslide.bounds-y property
+
+    Returns:
+        Feature with shifted polygon coordinates and updated centroid measurements.
+    """
+    if not feature or "geometry" not in feature:
+        return feature
+
+    geometry = feature["geometry"]
+    if "coordinates" not in geometry:
+        return feature
+
+    shifted_coordinates = []
+    for ring in geometry["coordinates"]:
+        shifted_ring = [
+            [point[0] - bounds_x, point[1] - bounds_y]
+            for point in ring
+        ]
+        shifted_coordinates.append(shifted_ring)
+    geometry["coordinates"] = shifted_coordinates
+
+    if "properties" in feature and "measurements" in feature["properties"]:
+        for measurement in feature["properties"]["measurements"]:
+            if measurement["name"] == "centroidX":
+                measurement["value"] -= bounds_x
+            elif measurement["name"] == "centroidY":
+                measurement["value"] -= bounds_y
+
+    return feature
+
+
 def deduplicate(features: list[dict], max_dist: float = 15 / 2) -> list[dict]:
     """
     Deduplicate cells based on their size and location.
@@ -705,7 +777,7 @@ def shapely_polygon_to_geojson(
     coords = np.array(polygon.exterior.coords.xy).T
     center = coords.mean(axis=1).tolist()
     coords = coords.tolist()
-    coords.append(coords[0])
+    coords.append(coords[0].copy())
     feature_id = str(uuid.uuid4())
     properties = {
         "objectType": object_type,
@@ -1241,6 +1313,15 @@ def main(args):
         logger.info("Filtering cells based on tissue contours")
         tissue_cnts = list(slide.tissue_cnts)
         polygons = filter_cells_by_contours(polygons, tissue_cnts)
+
+        bounds_x_val = float(slide.bounds_x.value)
+        bounds_y_val = float(slide.bounds_y.value)
+        if bounds_x_val != 0 or bounds_y_val != 0:
+            tissue_cnts = [
+                shapely.affinity.translate(cnt, xoff=-bounds_x_val, yoff=-bounds_y_val)
+                for cnt in tissue_cnts
+            ]
+
         tissue_features = []
         for i, cnt in enumerate(tissue_cnts):
             tissue_features.extend(
@@ -1311,6 +1392,15 @@ def main(args):
                     if polygon is not None:
                         artefact_polygons.append(polygon)
 
+
+        bounds_x_val = float(slide.bounds_x.value)
+        bounds_y_val = float(slide.bounds_y.value)
+        if bounds_x_val != 0 or bounds_y_val != 0:
+            artefact_polygons = [
+                shapely.affinity.translate(poly, xoff=-bounds_x_val, yoff=-bounds_y_val)
+                for poly in artefact_polygons
+            ]
+
         artefact_features = []
         for i, poly in enumerate(artefact_polygons):
             artefact_features.extend(
@@ -1340,6 +1430,15 @@ def main(args):
         with open(output_folder / artefact_contours_filename, "w") as f:
             json.dump(artefact_contours_fmt, f)
 
+    bounds_x_val = float(slide.bounds_x.value)
+    bounds_y_val = float(slide.bounds_y.value)
+    if bounds_x_val != 0 or bounds_y_val != 0:
+        logger.info(f"Applying bounds offset to output coordinates: x={bounds_x_val}, y={bounds_y_val}")
+        polygons = [
+            apply_bounds_offset_to_feature(f, bounds_x_val, bounds_y_val)
+            for f in polygons
+        ]
+
     geojson_fmt = {
         "type": "FeatureCollection",
         "features": polygons,

src/classpose/entrypoints/predict_wsi_cpsam.py

@@ -72,6 +72,7 @@
     DEFAULT_TRAIN_MPP,
     MIN_TILE_SIZE,
     SlideLoader,
+    apply_bounds_offset_to_feature,
     deduplicate,
     filter_cells_by_artefacts,
     filter_cells_by_contours,
@@ -434,6 +435,17 @@ def main(args):
         logger.info("Filtering cells based on tissue contours")
         tissue_cnts = list(slide.tissue_cnts)
         polygons = filter_cells_by_contours(polygons, tissue_cnts)
+
+        bounds_x_val = float(slide.bounds_x.value)
+        bounds_y_val = float(slide.bounds_y.value)
+        if bounds_x_val != 0 or bounds_y_val != 0:
+            tissue_cnts = [
+                shapely.affinity.translate(
+                    cnt, xoff=-bounds_x_val, yoff=-bounds_y_val
+                )
+                for cnt in tissue_cnts
+            ]
+
         tissue_features = []
         for i, cnt in enumerate(tissue_cnts):
             tissue_features.extend(
@@ -476,7 +488,7 @@ def main(args):
                 artefact_cnts,
                 artefact_geojson,
             ) = detect_artefacts_wsi(
-                slide=OpenSlide(args.slide_path),
+                slide=OpenSlide(slide.real_slide_path),
                 model_art_path=args.artefact_detection_model_path,
                 model_td_path=args.tissue_detection_model_path,
                 device=devices[0],
@@ -508,6 +520,16 @@ def main(args):
                     if polygon is not None:
                         artefact_polygons.append(polygon)
 
+        bounds_x_val = float(slide.bounds_x.value)
+        bounds_y_val = float(slide.bounds_y.value)
+        if bounds_x_val != 0 or bounds_y_val != 0:
+            artefact_polygons = [
+                shapely.affinity.translate(
+                    poly, xoff=-bounds_x_val, yoff=-bounds_y_val
+                )
+                for poly in artefact_polygons
+            ]
+
         artefact_features = []
         for i, poly in enumerate(artefact_polygons):
             artefact_features.extend(
@@ -539,6 +561,19 @@ def main(args):
         with open(output_folder / artefact_contours_filename, "w") as f:
             json.dump(artefact_contours_fmt, f)
 
+    bounds_x_val = float(slide.bounds_x.value)
+    bounds_y_val = float(slide.bounds_y.value)
+    if bounds_x_val != 0 or bounds_y_val != 0:
+        logger.info(
+            "Applying bounds offset to output coordinates: x=%s, y=%s",
+            bounds_x_val,
+            bounds_y_val,
+        )
+        polygons = [
+            apply_bounds_offset_to_feature(f, bounds_x_val, bounds_y_val)
+            for f in polygons
+        ]
+
     geojson_fmt = {
         "type": "FeatureCollection",
         "features": polygons,

src/classpose/grandqc/wsi_artefact_detection.py

@@ -62,6 +62,7 @@ def detect_artefacts_wsi(
     mpp_model_td: int = 10,
     m_p_s_model_td: int = 512,
     min_area: int = 0,
+    apply_bounds_offset: bool = False,
 ) -> tuple[np.ndarray, np.ndarray, dict[str, Any], dict[str, Any]]:
     """
     Detects artefacts in a whole-slide image using GrandQC method.
@@ -85,6 +86,9 @@ def detect_artefacts_wsi(
         mpp_model_td (int): MPP for tissue detection.
         m_p_s_model_td (int): patch size for tissue detection.
         min_area (int): minimum area for tissue polygons.
+        apply_bounds_offset (bool): if True, shift artefact contours and
+        GeoJSON by OpenSlide bounds offsets, so output coordinates are 
+        relative to the displayed image origin. Defaults to False.
 
     Returns:
         tuple: artefact_mask (np.ndarray), artefact_map (np.ndarray), artefact_cnts (dict), geojson (dict)
@@ -100,8 +104,12 @@ def detect_artefacts_wsi(
         encoder_model_weights=encoder_model_weights,
         device=device,
         min_area=min_area,
+        apply_bounds_offset=False,
     )
 
+    bounds_x = float(slide.properties.get("openslide.bounds-x", 0.0))
+    bounds_y = float(slide.properties.get("openslide.bounds-y", 0.0))
+
     # Load artefact model
     model_art_path = download_if_unavailable(model_art_path, MODEL_URL_PATH)
     preprocessing_fn = smp.encoders.get_preprocessing_fn(
@@ -306,6 +314,27 @@ def detect_artefacts_wsi(
         min_artifact_area,
     )
 
+    if apply_bounds_offset and (bounds_x != 0 or bounds_y != 0):
+        grandqc_logger.info(
+            "Applying bounds offset to artefact output coordinates: x=%s, y=%s",
+            bounds_x,
+            bounds_y,
+        )
+
+        for cnt in artefact_cnts.values():
+            cnt["contour"] = cnt["contour"] - np.array([bounds_x, bounds_y])
+            if cnt["holes"]:
+                cnt["holes"] = [
+                    hole - np.array([bounds_x, bounds_y]) for hole in cnt["holes"]
+                ]
+
+        for feature in geojson["features"]:
+            rings = feature["geometry"]["coordinates"]
+            feature["geometry"]["coordinates"] = [
+                [[point[0] - bounds_x, point[1] - bounds_y] for point in ring]
+                for ring in rings
+            ]
+
     del model
     del preprocessing_fn
     del slide
@@ -370,6 +399,7 @@ def detect_artefacts_wsi(
         model_td_path=args.model_td_path,
         min_area=args.min_area,
         device=args.device,
+        apply_bounds_offset=True,
     )
 
     # Save outputs

src/classpose/grandqc/wsi_tissue_detection.py

@@ -35,6 +35,7 @@ def detect_tissue_wsi(
     device: str = "cuda",
     # filtering parameters
     min_area: int = 0,
+    apply_bounds_offset: bool = False,
 ) -> tuple[Image, np.ndarray, np.ndarray, list[np.ndarray], dict[str, Any],]:
     """
     Detects tissue in a whole-slide image using a U-Net model. This is adapted
@@ -58,6 +59,9 @@ def detect_tissue_wsi(
         device (str, optional): device to run the model on. Defaults to "cuda".
         min_area (int, optional): minimum area of the tissue polygons. Defaults
             to 0.
+        apply_bounds_offset (bool, optional): if True, shift output contours and
+        GeoJSON by OpenSlide bounds offsets, so output coordinates are relative 
+        to the displayed image origin. Defaults to False.
 
     Returns:
         tuple: tuple containing the following elements:
@@ -86,6 +90,9 @@ def detect_tissue_wsi(
     model.to(device)
     model.eval()
 
+    bounds_x = float(slide.properties.get("openslide.bounds-x", 0.0))
+    bounds_y = float(slide.properties.get("openslide.bounds-y", 0.0))
+
     w_l0, h_l0, mpp, thumbnail_dimensions = extract_slide_info(
         slide, mpp_model_td
     )
@@ -265,6 +272,26 @@ def detect_tissue_wsi(
             }
         )
 
+    if apply_bounds_offset and (bounds_x != 0 or bounds_y != 0):
+        grandqc_logger.info(
+            "Applying bounds offset to tissue output coordinates: x=%s, y=%s",
+            bounds_x,
+            bounds_y,
+        )
+        for cnt in output_cnts.values():
+            cnt["contour"] = cnt["contour"] - np.array([bounds_x, bounds_y])
+            if cnt["holes"]:
+                cnt["holes"] = [
+                    hole - np.array([bounds_x, bounds_y]) for hole in cnt["holes"]
+                ]
+
+        for feature in geojson["features"]:
+            rings = feature["geometry"]["coordinates"]
+            feature["geometry"]["coordinates"] = [
+                [[point[0] - bounds_x, point[1] - bounds_y] for point in ring]
+                for ring in rings
+            ]
+
     del model
     del preprocessing_fn
     del slide
@@ -322,6 +349,7 @@ def detect_tissue_wsi(
         model_td_path=args.model_path,
         min_area=args.min_area,
         device=args.device,
+        apply_bounds_offset=True,
     )
 
     image.save(args.output_path + "_image.png")