diff --git a/app.py b/app.py index 6a59b0b..40fac35 100644 --- a/app.py +++ b/app.py @@ -14,6 +14,7 @@ import matplotlib as mpl import matplotlib.pyplot as plt from matplotlib import font_manager import numpy as np +from scipy import ndimage import streamlit as st from config import ( @@ -221,6 +222,52 @@ def plane_axis(plane: str) -> int: return {"Sagittal": 0, "Coronal": 1, "Axial": 2}[plane] +def plane_shift_from_xyz(plane: str, shift_xyz: Sequence[float]) -> Tuple[float, float]: + sx, sy, sz = (float(v) for v in shift_xyz) + if plane == "Axial": + return sx, sy + if plane == "Coronal": + return sx, sz + if plane == "Sagittal": + return sy, sz + raise ValueError(f"未知平面: {plane}") + + +def apply_manual_slice_transform( + image: np.ndarray, + plane: str, + shift_xyz: Sequence[float], + scale: float, + rotation_deg: float, +) -> np.ndarray: + """对被叠加切片做人工 2D 微调，仅用于可视化。""" + + scale = max(float(scale), 1e-3) + shift = np.asarray(plane_shift_from_xyz(plane, shift_xyz), dtype=np.float64) + theta = np.deg2rad(float(rotation_deg)) + forward = scale * np.asarray( + [ + [np.cos(theta), -np.sin(theta)], + [np.sin(theta), np.cos(theta)], + ], + dtype=np.float64, + ) + inverse = np.linalg.inv(forward) + center = (np.asarray(image.shape[:2], dtype=np.float64) - 1.0) / 2.0 + offset = center - inverse @ center - shift + return ndimage.affine_transform( + image.astype(np.float32, copy=False), + matrix=inverse, + offset=offset, + output_shape=image.shape, + order=1, + mode="constant", + cval=0.0, + prefilter=False, + output=np.float32, + ) + + def checkerboard(fixed_slice: np.ndarray, warped_slice: np.ndarray, tile: int = 24) -> np.ndarray: fixed_norm = normalize_slice(fixed_slice) warped_norm = normalize_slice(warped_slice) @@ -388,11 +435,13 @@ def render_three_way_views(fixed_xyz: np.ndarray, moving_xyz: np.ndarray, warped render_image(get_slice(volume, plane, index), f"{label} #{index}", plane=plane) -def render_overlay_views(fixed_xyz: np.ndarray, warped_xyz: np.ndarray) -> None: - fixed_xyz, warped_xyz = crop_to_common_shape(fixed_xyz, warped_xyz) +def render_overlay_views(fixed_xyz: np.ndarray, moving_xyz: np.ndarray, warped_xyz: np.ndarray) -> None: + fixed_xyz, moving_xyz, warped_xyz = crop_to_common_shape(fixed_xyz, moving_xyz, warped_xyz) + target_label = st.radio("叠加对象", ["配准后图像", "移动图像"], horizontal=True) + target_xyz = warped_xyz if target_label == "配准后图像" else moving_xyz mode = st.radio("对比模式", ["灰度融合", "红蓝差异", "差异热图", "棋盘格"], horizontal=True) if mode in {"灰度融合", "红蓝差异"}: - alpha = st.slider("配准后图像权重", 0.0, 1.0, 0.5, 0.05) + alpha = st.slider("被叠加图像权重", 0.0, 1.0, 0.5, 0.05) else: alpha = 0.5 if mode == "棋盘格": @@ -401,9 +450,41 @@ def render_overlay_views(fixed_xyz: np.ndarray, warped_xyz: np.ndarray) -> None: tile = 24 if mode == "红蓝差异": - st.caption("红/金色偏向配准后图像，蓝色偏向固定图像；白灰区域表示两者灰度接近。") + st.caption("红/金色偏向被叠加图像，蓝色偏向固定图像；白灰区域表示两者灰度接近。") elif mode == "差异热图": - st.caption("差异热图显示配准后图像减固定图像的灰度差，色条以当前切片的 99% 绝对差自动缩放。") + st.caption("差异热图显示被叠加图像减固定图像的灰度差，色条以当前切片的 99% 绝对差自动缩放。") + + with st.expander("人工微调重叠", expanded=True): + enabled = st.toggle("启用人工微调", value=False) + c1, c2, c3, c4 = st.columns(4) + with c1: + shift_x = st.slider("X 平移", -80.0, 80.0, 0.0, 0.5) + with c2: + shift_y = st.slider("Y 平移", -80.0, 80.0, 0.0, 0.5) + with c3: + shift_z = st.slider("Z 平移", -80.0, 80.0, 0.0, 0.5) + with c4: + scale = st.slider("缩放", 0.85, 1.15, 1.0, 0.005) + r1, r2, r3 = st.columns(3) + with r1: + axial_rotation = st.slider("轴状面旋转", -15.0, 15.0, 0.0, 0.5) + with r2: + coronal_rotation = st.slider("冠状面旋转", -15.0, 15.0, 0.0, 0.5) + with r3: + sagittal_rotation = st.slider("矢状面旋转", -15.0, 15.0, 0.0, 0.5) + if enabled: + st.caption( + f"当前人工微调：X={shift_x:+.1f}, Y={shift_y:+.1f}, Z={shift_z:+.1f}, " + f"缩放={scale:.3f}。这些参数只影响当前重叠显示。" + ) + + shift_xyz = (shift_x, shift_y, shift_z) if enabled else (0.0, 0.0, 0.0) + rotation_by_plane = { + "Axial": axial_rotation if enabled else 0.0, + "Coronal": coronal_rotation if enabled else 0.0, + "Sagittal": sagittal_rotation if enabled else 0.0, + } + scale_value = scale if enabled else 1.0 columns = st.columns(3) planes = [("Axial", "轴状面"), ("Coronal", "冠状面"), ("Sagittal", "矢状面")] @@ -412,15 +493,22 @@ def render_overlay_views(fixed_xyz: np.ndarray, warped_xyz: np.ndarray) -> None: with col: index = st.slider(f"{cn_name}切片", 0, fixed_xyz.shape[axis] - 1, fixed_xyz.shape[axis] // 2, key=f"overlay-{plane}") fixed_slice = get_slice(fixed_xyz, plane, index) - warped_slice = get_slice(warped_xyz, plane, index) + target_slice = get_slice(target_xyz, plane, index) + target_slice = apply_manual_slice_transform( + target_slice, + plane=plane, + shift_xyz=shift_xyz, + scale=scale_value, + rotation_deg=rotation_by_plane[plane], + ) if mode == "棋盘格": - render_image(checkerboard(fixed_slice, warped_slice, tile=tile), f"{cn_name} 棋盘格 #{index}", plane=plane) + render_image(checkerboard(fixed_slice, target_slice, tile=tile), f"{cn_name} 棋盘格 #{index}", plane=plane) elif mode == "红蓝差异": - render_rgb(color_difference_overlay(fixed_slice, warped_slice, alpha=alpha), f"{cn_name} 红蓝差异 #{index}", plane=plane) + render_rgb(color_difference_overlay(fixed_slice, target_slice, alpha=alpha), f"{cn_name} 红蓝差异 #{index}", plane=plane) elif mode == "差异热图": - render_signed_difference(signed_difference(fixed_slice, warped_slice), f"{cn_name} 差异热图 #{index}", plane=plane) + render_signed_difference(signed_difference(fixed_slice, target_slice), f"{cn_name} 差异热图 #{index}", plane=plane) else: - render_image(gray_fusion(fixed_slice, warped_slice, alpha=alpha), f"{cn_name} 灰度融合 #{index}", plane=plane) + render_image(gray_fusion(fixed_slice, target_slice, alpha=alpha), f"{cn_name} 灰度融合 #{index}", plane=plane) def render_ddf_views(fixed_xyz: np.ndarray, ddf_xyz: np.ndarray) -> None: @@ -830,7 +918,7 @@ def main() -> None: if warped_xyz is None: st.warning("尚未生成配准后图像。") else: - render_overlay_views(fixed_xyz, warped_xyz) + render_overlay_views(fixed_xyz, moving_xyz, warped_xyz) with tab_ddf: if ddf_xyz is None: