当前位置:首页>python>Python医学图像处理——光流法追踪特征点

Python医学图像处理——光流法追踪特征点

  • 2026-09-06 03:20:42
Python医学图像处理——光流法追踪特征点

在视频图像处理中,一个非常基础的问题是:当前帧中的某个结构,在下一帧移动到了哪里?

对于单张图像,我们可以做边缘检测、角点检测或目标分割,从而知道“这一帧中哪里有结构”。

但视频并不是一张张独立图像的简单堆叠。真正重要的是:

这些结构在相邻帧之间是否发生了移动?移动的方向是什么?移动了多远?

本质就是回答两件事情:

  1. 1. 这一帧中有什么特征?

  2. 2. 这些东西在下一帧去了哪里?

  3. 光流法(Optical Flow)正是为了解决这类问题而提出的。

光流的本质是建立在相邻帧之间局部图像结构的对应关系,并用一个运动向量描述每个位置的表观位移。

https://viso.ai/deep-learning/optical-flow/

光流是一种底层的运动估计方法。它可以作为物体追踪的基础:

当我们把属于同一个目标的多个局部运动向量组合起来时,就可以进一步推断这个目标的整体运动、局部形变或运动轨迹。但仅靠光流本身,通常还不能完整解决“这个物体是谁、是否仍然是同一个物体、是否发生遮挡或分裂”等更高层次的追踪问题。

https://opencv24-python-tutorials.readthedocs.io/en/latest/py_tutorials/py_video/py_lucas_kanade/py_lucas_kanade.html

光流法的具体实现可以分为两个步骤:

  1. 1. 先在当前帧中选出一些容易被稳定追踪的点,例如角点或纹理丰富的区域,

  2. 2. 然后利用光流法估计这些点在下一帧中的新位置。


  1. 这篇文章会一段冠脉造影数据作为例子,通过 Python 实现Shi-Tomasi 特征点检测 + Lucas-Kanade 光流追踪。具体效果:

数据来源:https://www.youtube.com/shorts/4pkwqjj6ZSg

一、原理解析

1、什么是光流法

光流法(Optical Flow)用于估计连续图像帧之间的运动。

更直观地说,假设在第 t 帧中有一个局部图像结构:

第 t 帧:      *

在第 t+1 帧中,它移动到了右下方:

第 t+1 帧:        *

那么光流法希望估计的就是这个运动向量:

也就是:这个点从哪里来,又移动到了哪里。

根据估计对象的不同,光流一般可以分成两类:

  1. 1. 稠密光流(Dense Optical Flow):对图像中几乎每个像素都估计运动。

  2. 2. 稀疏光流(Sparse Optical Flow):只选取一部分可靠的特征点,然后追踪这些点。

本文使用的是第二种:稀疏光流。

因为在很多科研图像中,我们并不一定需要知道每个像素的运动,只需要追踪一些稳定、可靠、容易识别的局部点,就可以很好地观察整体运动趋势。


2、光流法的基本原理

光流法最核心的假设是:同一个局部图像结构,在很短时间内亮度近似不变,只是位置发生了变化。

用公式表示就是:

其中:

这就是所谓的亮度恒定假设(Brightness Constancy)。

但是,只靠一个像素是不够的。

因为一个像素只有一个灰度值,而我们要估计的是两个未知量:

所以 Lucas-Kanade 光流法的思路是:

不只看一个像素,而是看这个像素周围的一小块窗口,并假设这个窗口里的局部结构整体一起移动。

例如,对于一个特征点,不是只比较这个点本身,而是比较它周围 21 x 21 的小区域在下一帧中移动到了哪里。

这就把问题从“一个像素的灰度变化”变成了“一个小图像块的匹配问题”。

在 OpenCV 中,Lucas-Kanade 光流对应的函数是:

cv2.calcOpticalFlowPyrLK(...)

其中 PyrLK 中的 Pyr 表示 Pyramid,也就是图像金字塔。

为什么需要金字塔?

因为相邻两帧之间的运动可能并不总是很小。如果直接在原始分辨率上追踪,点可能会跳出局部搜索范围。图像金字塔会先在低分辨率图像上估计粗略位移,再回到高分辨率图像上细化结果。

可以理解为:

低分辨率:先找大概移动方向高分辨率:再修正具体位置

3、特征点的选取原理

在使用稀疏光流时,一个关键问题是: 应该追踪哪些点?

并不是所有像素都适合被追踪。

1、平坦区域不适合追踪

如果一个区域几乎没有纹理,例如一大片灰色背景:

灰 灰 灰 灰灰 灰 灰 灰灰 灰 灰 灰

那么其中任意一个点移动后,看起来都差不多。算法很难判断它到底移动到了哪里。

2、单纯边缘也不够稳定

如果一个点位于一条边缘上:

暗 暗 | 亮 亮暗 暗 | 亮 亮暗 暗 | 亮 亮

它在垂直于边缘的方向上比较容易定位,但沿着边缘方向移动时会产生歧义。

这就是经典的孔径问题(Aperture Problem)。

3、角点更适合追踪

角点附近在两个方向上都有明显灰度变化:

暗 暗 | 亮 亮暗 暗 | 亮 亮------+------亮 亮 | 暗 暗亮 亮 | 暗 暗

这类点局部结构更有辨识度,因此更适合作为追踪点。

Shi-Tomasi 角点检测就是一种经典的特征点选取方法。它的基本思想是:考察一个局部窗口在不同方向上的灰度变化。如果这个窗口在两个方向上都有足够明显的变化,那么它就是一个好的特征点。

https://www.datatechnotes.com/2023/10/shi-tomasi-corner-detection-example.html

在 OpenCV 中,对应函数是:

cv2.goodFeaturesToTrack(...)

本文使用 Shi-Tomasi 检测关键帧中的特征点,然后使用 Lucas-Kanade 光流追踪这些点。

整体流程如下:

关键帧  ↓Shi-Tomasi 角点检测  ↓Lucas-Kanade 光流逐帧追踪  ↓过滤错误匹配和异常大位移  ↓绘制运动轨迹

二、Python 环境安装

本文只需要两个主要 Python 包:

opencv-pythonnumpy

可以直接安装:

pip install opencv-python numpy

也可以新建一个 conda 环境:

conda create -n opticalflow python=3.10conda activate opticalflowpip install opencv-python numpy

导入需要的包:

from pathlib import Pathimport cv2import numpy as np

三、实例代码

1、读取图像

from pathlib import Pathimport cv2import numpy as npbase_dir = Path(__file__).resolve().parentframe_dir = base_dir / "Coronary_angiography_video_tif_8bit"out_dir = base_dir / "optical_flow_tutorial_result"each_frame_dir = base_dir / "optical_flow_tracking_frames_colored"out_dir.mkdir(exist_ok=True)each_frame_dir.mkdir(exist_ok=True)start_frame = 1end_frame = 11def read_gray(path):    image = cv2.imdecode(np.fromfile(str(path), dtype=np.uint8), cv2.IMREAD_GRAYSCALE)    if image is None:        raise FileNotFoundError(path)    return imagedef save_image(path, image):    ok, buffer = cv2.imencode(path.suffix, image)    if not ok:        raise RuntimeError(f"Cannot save image: {path}")    buffer.tofile(str(path))def frame_path(index):    return frame_dir / f"frame_{index:06d}.tif"prev_gray = read_gray(frame_path(start_frame))save_image(out_dir / "01_start_frame.png", prev_gray)

读取到的第 1 帧就是前面展示的起始图像。

2、Shi-Tomasi 检测特征点

接下来在关键帧上检测适合追踪的角点:

feature_params = dict(    maxCorners=180,    qualityLevel=0.01,    minDistance=7,    blockSize=7,)points = cv2.goodFeaturesToTrack(prev_gray, mask=None, **feature_params)

这些参数的含义是:

  • maxCorners=180:最多检测 180 个点;

  • qualityLevel=0.01:角点质量阈值;

  • minDistance=7:两个特征点之间的最小距离;

  • blockSize=7:计算角点响应时使用的局部窗口大小。

把检测到的点画出来:

def make_random_colormap(count, seed=19):    rng = np.random.default_rng(seed)    hues = np.linspace(0, 179, count, endpoint=False, dtype=np.float32)    rng.shuffle(hues)    hsv = np.zeros((count, 1, 3), dtype=np.uint8)    hsv[:, 0, 0] = hues.astype(np.uint8)    hsv[:, 0, 1] = rng.integers(170, 245, size=count, dtype=np.uint8)    hsv[:, 0, 2] = rng.integers(210, 255, size=count, dtype=np.uint8)    return cv2.cvtColor(hsv, cv2.COLOR_HSV2BGR).reshape(count, 3)def draw_points(gray, points, colors, point_ids):    canvas = cv2.cvtColor(gray, cv2.COLOR_GRAY2BGR)    for point, point_id in zip(points.reshape(-1, 2), point_ids):        x, y = point        color = colors[int(point_id)].tolist()        cv2.circle(canvas, (int(x), int(y)), 3, color, -1, cv2.LINE_AA)        cv2.circle(canvas, (int(x), int(y)), 4, (255, 255, 255), 1, cv2.LINE_AA)    return canvastrack_colors = make_random_colormap(len(points))initial_ids = np.arange(len(points))save_image(out_dir / "02_shi_tomasi_points.png", draw_points(prev_gray, points, track_colors, initial_ids))save_image(each_frame_dir / "tracked_000040.png", draw_points(prev_gray, points, track_colors, initial_ids))

结果如下:

可以看到,特征点主要分布在血管边缘、分叉附近以及局部灰度变化明显的位置。

这也符合前面的原理:光流追踪更喜欢局部结构明显的点。

3、Lucas-Kanade 光流追踪

定义 Lucas-Kanade 光流参数:

lk_params = dict(    winSize=(21, 21),    maxLevel=3,    criteria=(cv2.TERM_CRITERIA_EPS | cv2.TERM_CRITERIA_COUNT, 30, 0.01),)

其中:

  • winSize=(21, 21):使用 21 x 21 的局部窗口进行匹配;

  • maxLevel=3:使用 3 层图像金字塔;

  • criteria:迭代停止条件。

对相邻两帧进行光流估计:

curr_gray = read_gray(frame_path(41))next_points, status, error = cv2.calcOpticalFlowPyrLK(    prev_gray,    curr_gray,    points,    None,    **lk_params,)

这里:

  • points 是关键帧中的特征点;

  • next_points 是这些点在下一 帧中的位置;

  • status 表示每个点是否追踪成功;

  • error 表示局部窗口匹配误差。

位移向量如下:

每一个箭头表示一个特征点从第1帧到第2帧的移动方向。

4、过滤错误匹配和异常大位移

实际追踪时,不能只看 status。

有些点虽然被 OpenCV 返回为追踪成功,但可能跳到了错误位置,表现为一个特别长的箭头。

因此这里加入四个过滤条件:

status = status.reshape(-1).astype(bool)error = error.reshape(-1)displacement = np.linalg.norm(    next_points.reshape(-1, 2) - points.reshape(-1, 2),    axis=1,)h, w = curr_gray.shapexy = next_points.reshape(-1, 2)in_bounds = (    (xy[:, 0] >= 0) & (xy[:, 0] < w) &    (xy[:, 1] >= 0) & (xy[:, 1] < h))good = (    status &    (error < 25.0) &    (displacement <= 30.0) &    in_bounds)

这里的过滤逻辑是:

  • status:OpenCV 认为追踪成功;

  • error < 25.0:局部窗口匹配误差不能太大;

  • displacement <= 30.0:单帧位移不能异常大;

  • in_bounds:点仍然在图像范围内。

加入单帧位移阈值后,可以去掉少数明显不合理的跳变点,避免最终图中出现特别长的错误箭头。

5、逐帧追踪并绘制轨迹

从第 1 帧追踪到第 10 帧,可以写成一个循环:

tracks = [[tuple(p.ravel())] for p in points]active_ids = np.arange(len(points))max_displacement = 30.0rejected_by_displacement = 0for frame_index in range(start_frame + 1, end_frame + 1):    curr_gray = read_gray(frame_path(frame_index))    next_points, status, error = cv2.calcOpticalFlowPyrLK(        prev_gray,        curr_gray,        points,        None,        **lk_params,    )    status = status.reshape(-1).astype(bool)    error = error.reshape(-1)    displacement = np.linalg.norm(next_points.reshape(-1, 2) - points.reshape(-1, 2), axis=1)    h, w = curr_gray.shape    xy = next_points.reshape(-1, 2)    in_bounds = (xy[:, 0] >= 0) & (xy[:, 0] < w) & (xy[:, 1] >= 0) & (xy[:, 1] < h)    good = status & (error < 25.0) & (displacement <= max_displacement) & in_bounds    rejected_by_displacement += int(np.count_nonzero(        status & (error < 25.0) & in_bounds & (displacement > max_displacement)    ))    good_new = next_points[good]    good_ids = active_ids[good]    for track_id, point in zip(good_ids, good_new):        tracks[int(track_id)].append(tuple(point.ravel()))    points = good_new.reshape(-1, 1, 2).astype(np.float32)    active_ids = good_ids    prev_gray = curr_gray

最后绘制轨迹:

def draw_tracks(gray, tracks, colors):    canvas = cv2.cvtColor(gray, cv2.COLOR_GRAY2BGR)    for i, track in enumerate(tracks):        if len(track) < 2:            continue        pts = np.array(track, dtype=np.int32)        color = colors[i].tolist()        cv2.polylines(canvas, [pts], False, color, 1, cv2.LINE_AA)        cv2.circle(canvas, tuple(pts[-1]), 2, color, -1, cv2.LINE_AA)    return canvasfinal_gray = read_gray(frame_path(end_frame))save_image(out_dir / "04_tracks_40_to_50.png", draw_tracks(final_gray, tracks, track_colors))

最终追踪结果如下:

6、结果分析

运行脚本后,输出结果为:

Initial feature points: 180Full tracks from frame 1 to 11: 158Rejected large jumps: 4

也就是说:

  • 第 1 帧初始检测到 180 个 Shi-Tomasi 特征点;

  • 从第 1 帧追踪到第 11 帧后,有 158 条轨迹完整保留下来;

  • 有 4 次追踪结果因为单帧位移过大,被认为是不合理跳变并剔除。

从最终轨迹图中可以看到,很多轨迹集中在血管边缘和分叉区域附近。 这是因为这些区域灰度变化明显,更容易被 Shi-Tomasi 检测为特征点,也更容易被 Lucas-Kanade 光流追踪。

本篇教程的所有Python代码,可以在GitHub上下载:

https://github.com/ethanzhao9/Medical-Image-Processing

最新文章

随机文章