如何检测将显微镜图像旋转对齐到模板的良好特征

Posted 2023-03-16

【中文标题】如何检测将显微镜图像旋转对齐到模板的良好特征

【英文标题】：How to detect good features for rotationally aligning microscope images to a template

【发布时间】：2019-11-27 16:23:32

【问题描述】：

我正在开展一个项目，以自动旋转流体实验的显微镜图像堆栈，以便它们与微流控芯片的 CAD 模板图像对齐。我正在使用 Python 中的 OpenCV 包进行图像处理。必须有正确的旋转方向，以便可以正确屏蔽图像以进行分析。我们的芯片具有填充有荧光染料的标记，在每一帧中都可见。模板和示例图像如下所示（模板可以缩放到任意大小，但图像的相关区域通常约为 100x100 像素左右）：

我无法将图像旋转对齐到 CAD 模板。通常，CAD 模板和图像之间的错位小于几度，这仍然足以干扰分析，所以我需要能够测量旋转差异，即使它相对较小。

以下在线示例我正在使用以下程序：

使用三次插值将图像放大到与模板大致相同的大小 (~800 x 800) 使用 Otsu 方法对两个图像进行阈值处理 使用内置方法（我尝试过 ORB、AKAZE 和 Brief）查找关键点并提取描述符。 使用具有汉明距离的蛮力匹配器匹配描述符。 获取最佳匹配并使用它们计算部分仿射变换矩阵 使用该矩阵推断旋转偏移，将一个图像扭曲到另一个图像作为检查。

这是我的代码示例（部分借用自 here）：

import numpy as np
import cv2
import matplotlib.pyplot as plt

MAX_FEATURES = 500
GOOD_MATCH_PERCENT = 0.5

def alignImages(im1, im2,returnpoints=False):
    # Detect ORB features and compute descriptors.
    size1 = int(0.1*(np.mean(np.shape(im1))))
    size2 = int(0.1*(np.mean(np.shape(im2))))
    orb1 = cv2.ORB_create(MAX_FEATURES,edgeThreshold=size1,patchSize=size1)
    orb2 = cv2.ORB_create(MAX_FEATURES,edgeThreshold=size2,patchSize=size2)
    keypoints1, descriptors1 = orb1.detectAndCompute(im1, None)
    keypoints2, descriptors2 = orb2.detectAndCompute(im2, None)

    matcher = cv2.BFMatcher(cv2.NORM_HAMMING,crossCheck=True)
    matches = matcher.match(descriptors1,descriptors2)
    # Sort matches by score
    matches.sort(key=lambda x: x.distance, reverse=False)

    # Remove not so good matches
    numGoodMatches = int(len(matches) * GOOD_MATCH_PERCENT)
    matches = matches[:numGoodMatches]

    # Draw top matches
    imMatches = cv2.drawMatches(im1, keypoints1, im2, keypoints2, matches, None)
    cv2.imwrite("matches.jpg", imMatches)

    # Extract location of good matches
    points1 = np.zeros((len(matches), 2), dtype=np.float32)
    points2 = np.zeros((len(matches), 2), dtype=np.float32)
    for i, match in enumerate(matches):
        points1[i, :] = keypoints1[match.queryIdx].pt
        points2[i, :] = keypoints2[match.trainIdx].pt
    # Find homography
    M, inliers = cv2.estimateAffinePartial2D(points1,points2)
    height, width = im2.shape
    im1Reg = cv2.warpAffine(im1,M,(width,height))
    return im1Reg, M

if __name__ == "__main__":
    test_template = cv2.cvtColor(cv2.imread("test_CAD_cropped.png"),cv2.COLOR_RGB2GRAY)
    test_image = cv2.cvtColor(cv2.imread("test_CAD_cropped.png"),cv2.COLOR_RGB2GRAY)
    fx = fy = 88/923
    test_image_big = cv2.resize(test_image,(0,0),fx=1/fx,fy=1/fy,interpolation=cv2.INTER_CUBIC)
    ret, imRef_t = cv2.threshold(test_template,0,255,cv2.THRESH_BINARY+cv2.THRESH_OTSU)
    ret, test_big_t = cv2.threshold(test_image_big,0,255,cv2.THRESH_BINARY+cv2.THRESH_OTSU)
    imReg, M = alignImages(test_big_t,imRef_t)
    fig, ax = plt.subplots(nrows=2,ncols=2,figsize=(8,8))
    ax[1,0].imshow(imReg)
    ax[1,0].set_title("Warped Image")
    ax[0,0].imshow(imRef_t)
    ax[0,0].set_title("Template")
    ax[0,1].imshow(test_big_t)
    ax[0,1].set_title("Thresholded Image")
    ax[1,1].imshow(imRef_t - imReg)
    ax[1,1].set_title("Diff")
    plt.show()

在这个例子中，我得到了以下错误的转换，因为只有 3 个匹配的关键点并且它们都是不正确的：

我发现无论我的关键点/描述符参数如何，我都倾向于获得太少的“好”功能。 我可以做些什么来更好地预处理我的图像以更可靠地获得好的特征，或者有没有更好的方法来将我的图像与这个不涉及关键点匹配的模板对齐？具体这个实验的应用意味着我不能使用获得专利的关键点提取器/描述符，如 SURF 和 SIFT。

【问题讨论】：

你试过Hough Transform来检查旋转差异吗？

【参考方案1】：

Fourier Mellin transform 是一种仅基于旋转、平移和缩放来对齐两个图像的好方法。

这是在DIPlib 中使用implementation 的示例（披露：我是作者）：

import diplib as dip

# load data
image = dip.ImageRead('image.png')
template = dip.ImageRead('template.png')
template = template.TensorElement(0) # this one is RGB, take any one channel

# pad the two images with zeros so they have equal sizes
sz = [max(image.Size(0), template.Size(0)), max(image.Size(1), template.Size(1))]
image = image.Pad(sz)
template = template.Pad(sz)

# match
res = dip.FourierMellinMatch2D(template, image)

# display
dip.JoinChannels((template,res,res)).Show()

---

但是，还有许多其他方法。这里的关键是模板和图像都非常简单，并且非常相似。这使得注册非常容易。

例如，假设您对模板进行了适当的缩放（我认为这应该不是问题），您需要做的就是找到旋转和平移。您可以强制旋转，只需将图像旋转一组小角度，并将每个结果与模板匹配（互相关）。具有最佳匹配（最大互相关值）的那个具有适当的旋转。如果您需要进行非常精确的旋转估计，您可以进行第二组角度接近第一组中的最佳选择。

互相关便宜且易于计算，并且可以产生高精度的平移估计（傅里叶梅林方法广泛使用它）。不要只在互相关输出中找到具有最大值的像素，您可以将抛物线拟合到该像素周围的几个像素上，并使用拟合抛物线的最大值的位置。这会导致平移的亚像素估计。