《Python计算机视觉》——局部图像描述子

2.1 Harris角点检测器

角点检测

Harris角点检测算法（也称为Harris&Stephens角点检测器）是一个极为简单的角点检测算法。该算法的主要思想是：如果像素周围显示存在多余一个方向的边，我们认为该点为兴趣点。该点称为角点。

我们把图像领域中的点x上的对称半正定矩阵$M_I=M_I(x)$定义为：

其中$\nabla I$为包含导数$I_x$和$I_y$的图像梯度。由于该定义，$M_I$的秩为1，特征值为${\lambda }_1=|\nabla I{|}^2$和${\lambda }_2=0$。现在对于图像的每一个像素，我们可以计算出该矩阵。

选择权重矩阵$W$（通常为高斯滤波器$G_{\sigma }$），我们可以得到卷积：

该卷积的目的是得到$M_I$在周围像素上的局部平均。计算出的矩阵$\overline{M_I}$有称为Harris矩阵。$M$的宽度决定了在像素x周围的感兴趣区域。像这样的区域附近对矩阵$\overline{M_I}$取平均的原因是，特征值会依赖于局部图像特性而变化。如果图像的梯度在该区域变化，那么$\overline{M_I}$的第二个特征值将不再为0。如果图像的梯度没有变化，$\overline{M_I}$的特征值也不会变化。
取决于该区域$\nabla I$的值，Harris矩阵$\overline{M_I}$的特征值有三种情况：

• 如果${\lambda }_1$和${\lambda }_2$都是很大的正数，则该x点为角点；
• 如果${\lambda }_1$很大。${\lambda }_2\approx 0$，则该区域内存在一个边，该区域内的平均$\overline{M_I}$的特征值不会变化太大；
• 如果${\lambda }_1\approx {\lambda }_2\approx 0$，该区域内为空。
  在不需要实际计算特征值的情况下，为了把重要的情况和其他情况分开，Harris和Stephens在文献[12]中引入了指示函数：
  $$\det (\overline{M_I})-\kappa trace(\overline{M_I}{)}^2$$
  为了去除加权常数$\kappa$，我们通常使用商数：
  $$\frac{det(\overline{M_I})}{trace(\overline{M_I}{)}^2}$$
  作为指示器。
  下面我们写出Harris角点检测程序：具体实现分为三个函数。

首先定义compute_harris_response() 函数，用于计算输入图像中每个像素点的角点响应值。该函数使用高斯导数实现，参数$\sigma$定义了使用的高斯滤波器的尺度大小。你也可以修改这个函数，对$x$和$y$方向上的尺度参数，以及尝试平均操作中的不同尺度，来计算Harris矩阵：

# # 角点检测
def compute_harris_response(im, sigma=3):
    """ 在一幅灰度图像中，对每个像素计算Harris角点检测器响应函数 """

    # 计算导数
    imx = zeros(im.shape)     # 创建一个与输入图像（im）相同大小的全零矩阵（imx）
    filters.gaussian_filter(im, (sigma, sigma), (0, 1), imx)   # 将输入图像（im）应用高斯滤波器，生成一个输出矩阵（imx），(0, 1)指定了滤波器的方向（沿图像的水平方向进行滤波）
    imy = zeros(im.shape)
    filters.gaussian_filter(im, (sigma, sigma), (1, 0), imy)   # (0, 1)滤波器的方向为垂直方向
    
    # 计算Harris矩阵的分量
    Wxx = filters.gaussian_filter(imx * imx, sigma)
    Wxy = filters.gaussian_filter(imx * imy, sigma)
    Wyy = filters.gaussian_filter(imy * imy, sigma)
    
    # 计算特征值和迹
    Wdet = Wxx * Wyy - Wxy ** 2     # 计算图像中每个位置的Hessian矩阵行列式值
    Wtr = Wxx + Wyy    # 计算每个位置上Hessian矩阵的迹，检测出图像中的边缘
    
    return Wdet / Wtr   # 返回图像中每个位置的角点响应值

上面的函数会返回像素值Harris响应函数值的一幅图像。现在，我们需要从这幅图像中挑选出需要的信息。然后，选取像素值高于阈值的所有图像点；再加上额外的限制，即角点之间的间隔必须大于设定的最小距离。这种方法会产生很好的角点检测结果。为了实现该算法，我们获取所有的候选像素点，以角点响应值递减的顺序排序，然后将距离已标记为角点位置过近的区域从候选像素点中删除。

其次，定义get_harris_points() 函数用于从 Harris 角点响应矩阵中提取角点。

def get_harris_points(harrisim, min_dist=10, threshold=0.1):
    """ 从一幅Harris响应图像中返回角点。min_dist为分割角点和图像边界的最小像素数目 """

    corner_threshold = harrisim.max() * threshold      # 用于将整个角点响应矩阵中的角点分为两类的阈值
    harrisim_t = (harrisim > corner_threshold) * 1     # 大于阈值的元素设置为1，其余元素设置为0，存放在矩阵harrisim_t中
    coords = array(harrisim_t.nonzero()).T             # 获取harrisim_t矩阵中所有非零元素（即角点）的坐标并转置
    candidate_values = [harrisim[c[0], c[1]] for c in coords]     # 获取harrisim_t矩阵中所有角点的响应值，并保存到一个列表中
    index = argsort(candidate_values)[::-1]      # 从大到小进行排序
    allowed_locations = zeros(harrisim.shape)    # 创建一个与harrisim形状相同的全零矩阵
    allowed_locations[min_dist:-min_dist, min_dist:-min_dist] = 1
    filtered_coords = []
    for i in index:
        if allowed_locations[coords[i, 0], coords[i, 1]] == 1:
            filtered_coords.append(coords[i])
            allowed_locations[(coords[i, 0] - min_dist):(coords[i, 0] + min_dist),
            (coords[i, 1] - min_dist):(coords[i, 1] + min_dist)] = 0
    return filtered_coords

现在你有了检测图像中角点所需要的所有函数。为了显示图像中的角点，你可以使用Matplotlib模块绘制函数。主要代码如下：

def plot_harris_points(image, filtered_coords):
    """ 绘制图像中检测到的角点"""

    figure()
    gray()
    imshow(image)
    plot([p[1] for p in filtered_coords],
         [p[0] for p in filtered_coords], '*')
    axis('off')
    show()

下面只需要在代码中引入需要用到的包，并写下以下代码就可以实现Harris角点检测。

from pylab import *
from numpy import *
from scipy.ndimage import filters
from PIL import Image

im = array(Image.open(r"course/computer_vision/data/all_data/1.pic.jpg").convert('L'))
harrisim = compute_harris_response(im)
filtered_coords = get_harris_points(harrisim)
plot_harris_points(im, filtered_coords)

整个过程为：首先，打开该图像，转换为灰度图像。然后，计算响应函数，基于响应值选择角点。最后，在原始图像中覆盖绘制检测出的角点。

输入原始图片为：

运行结果为：

分析：需要注意的是代码部分中，照片的路径建议选择绝对路径，避免运行出错。使用绝对路径比使用相对路径更可靠，更不容易出错。

在图像间寻找对应点

Harris角点检测器仅仅能够检测出图像中的兴趣点，但是没有给出通过比较图像间的兴趣点来寻找匹配角点的方法。我们需要在每个点上加入描述子信息，并给出一个比较这些描述子的方法。

兴趣点描述子是分配给兴趣点的一个向量，描述该点附近的图像的表观信息。描述子越好，寻找到的对应点越好。我们用对应点或者点的对应来描述相同物体和场景点在不同图像上形成的像素点。

Harris角点的描述子通常是由周围图像像素块的灰度值，以及用于比较的归一化互相关矩阵构成的。图像的像素块由以该像素点为中心的周围矩形部分图像构成。

通常，两个（相同大小）像素块$I_1(x)$和$I_2(x)$的相关矩阵定义为：
$$c(I_1,I_2)=\sum _{x}f(I_1(x),I_2(x))$$
其中，函数$f$随着相关方法的变化而变化。上式取像素位置$x$的和。对于互相关矩阵，函数$f(I_1,I_2)=I_1{I}_2$，因此，$c(I_1,I_2)=I_1\cdot I_2$，其中·表示向量乘法（按照行或者列堆积的像素）。$c(I_1,I_2)$的值越大，像素块$I_1$和$I_2$的相似度越高。

归一化的互相关矩阵是互相关矩阵的一种变形，可以定义为：
$$ncc(I_1,I_2)=\frac{1}{n-1}\sum _x\frac{(I_1(x)-{\mu }_1)}{{\sigma }_1}\frac{(I_2(x)-{\mu }_2)}{{\sigma }_2}$$
其中，$n$为像素块中像素的数目，${\mu }_1$和${\mu }_2$表示每个像素块中的平均像素值强度，${\sigma }_1$和${\sigma }_2$分别表示每个像素块中的标准差。通过减去均值和除以标准差，该方法对图像亮度变化具有稳健性。

为获取图像像素块，并使用归一化的互相关矩阵来比较它们，你需要另外两个函数，主要代码为下：

def get_descriptors(image, filtered_coords, wid=5):
    """ 对于每个返回的点，返回点周围2*wid+1个像素的值（假设选取的点min_distance>wid） """

    desc = []
    for coords in filtered_coords:
        patch = image[coords[0] - wid:coords[0] + wid + 1,
                coords[1] - wid:coords[1] + wid + 1].flatten()
        desc.append(patch)

    return desc

def match(desc1, desc2, threshold=0.5):
    """ 对于第一幅图像中的每个角点描述子，使用归一化互相关，选取它在第二幅图像中的匹配角点 """

    n = len(desc1[0])

    # 点对的距离
    d = -ones((len(desc1), len(desc2)))
    for i in range(len(desc1)):
        for j in range(len(desc2)):
            d1 = (desc1[i] - mean(desc1[i])) / std(desc1[i])
            d2 = (desc2[j] - mean(desc2[j])) / std(desc2[j])
            ncc_value = sum(d1 * d2) / (n - 1)
            if ncc_value > threshold:
                d[i, j] = ncc_value

    ndx = argsort(-d)
    matchscores = ndx[:, 0]

    return matchscores

第一个函数的参数为奇数大小长度的方形灰度图像块，该图像块的中心为处理的像素点。该函数将图像块像素值压平成一个向量，然后添加到描述子列表中。第二个函数使用归一化的互相关矩阵，将每个描述子匹配到另一个图像中的最优的候选点。由于数值较高的距离代表两个点能够更好的匹配，所以在排序之前，我们对距离取相反数。为了获得更稳定的匹配，我们从第二幅图像向第一幅图像匹配，然后过滤掉在两种方法中不都是最好的匹配。下面的函数可以实现该操作：

def match_twosided(desc1, desc2, threshold=0.5):
    """ 两边对称版本的match() """

    matches_12 = match(desc1, desc2, threshold)
    matches_21 = match(desc2, desc1, threshold)

    ndx_12 = where(matches_12 >= 0)[0]

    # 去除非对称的匹配
    for n in ndx_12:
        if matches_21[matches_12[n]] != n:
            matches_12[n] = -1

    return matches_12

这些匹配可以通过在两边分别绘制出图像，使用线段连接匹配的像素点来直观地可视化。下面的代码可以实现匹配点的可视化：

def appendimages(im1, im2):
    """ 返回将两幅图像并排拼接成的一幅新图像 """

    # 选取具有最少行数的图像，然后填充足够的空行
    rows1 = im1.shape[0]
    rows2 = im2.shape[0]

    if rows1 < rows2:
        im1 = concatenate((im1, zeros((rows2 - rows1, im1.shape[1]))), axis=0)
    elif rows1 > rows2:
        im2 = concatenate((im2, zeros((rows1 - rows2, im2.shape[1]))), axis=0)
    # 如果这些情况都没有，那么他们的行数相同，不需要填充

    return concatenate((im1, im2), axis=1)


def plot_matches(im1, im2, locs1, locs2, matchscores, show_below=True):
    """ 显示一幅带有连接匹配之间连线的图片
        输入：im1，im2（数组图像），locs1，locs2（特征位置），matchscores（match()的输出），
        show_below（如果图像应该显示在匹配的下方） """

    im3 = appendimages(im1, im2)
    if show_below:
        im3 = vstack((im3, im3))

    imshow(im3)

    cols1 = im1.shape[1]
    for i, m in enumerate(matchscores):
        if m > 0:
            plot([locs1[i][1], locs2[m][1] + cols1], [locs1[i][0], locs2[m][0]], 'c')
    axis('off')

为使用归一化的互相关矩阵来寻找对应点的例子，该图像可以通过下面的命令实现：

im1 = array(Image.open(r"course/computer_vision/data/all_data/5.pic.jpg").convert("L"))
im2 = array(Image.open(r"course/computer_vision/data/all_data/6.pic.jpg").convert("L"))

wid = 5
harrisim = compute_harris_response(im1, 5)
filtered_coords1 = get_harris_points(harrisim, wid + 1)
d1 = get_descriptors(im1, filtered_coords1, wid)

harrisim = compute_harris_response(im2, 5)
filtered_coords2 = get_harris_points(harrisim, wid + 1)
d2 = get_descriptors(im2, filtered_coords2, wid)

print('starting mathching')
matches = match_twosided(d1, d2)

figure()
gray()
plot_matches(im1, im2, filtered_coords1, filtered_coords2, matches)
show()

以上代码的运行结果为：

为了看的更清楚，你可以画出匹配的子集。在上面的代码中，可以通过将数组matches替换成matches[:100]或者任意子集来实现。

2.2 SIFT（尺度不变特征变换）

SIFT特征包括兴趣点检测器和描述子。SIFT描述子具有非常强的稳健性，这在很大程度上也是SIFT特征能够成功和流行的主要原因。SIFT特征对于尺度、旋转和亮度都具有不变性，因此，它可以用于三维视角和噪声的可靠匹配。

2.2.1 兴趣点

SIFT特征使用高斯差分函数来定位兴趣点：
$$D(x,\sigma )=[G_{\kappa \sigma }(x)-G_{\sigma }(x)]\ast I(x)=[G_{\kappa \sigma }-G_{\sigma }]\ast I=I_{\kappa \sigma }-I_{\sigma }$$
其中，$G_{\sigma }$是二维高斯核，$I_{\sigma }$是使用$G_{\sigma }$模糊的灰度图像，$\kappa$是决定相差尺度的常数。兴趣点是在图像位置和尺度变化下$D(x,\sigma )$的最大值和最小值点。这些候选位置点通过滤波去除不稳定。基于一些准则，比如认为低对比度和位于边上的点不是兴趣点，我们可以去除一些候选兴趣点。

2.2.2 描述子

上面讨论的兴趣点（关键点）位置描述子给出了兴趣点的位置和尺度信息。为了实现旋转不变性，基于每个点周围图像梯度的方向和大小，SIFT描述子又引入参考方向。SIFT描述子使用主方向描述参考方向。主方向使用方向直方图（以大小为权重）来度量。

下面我们基于位置、尺度和方向信息来计算描述子。为了对图像亮度具有稳健性，SIFT描述子使用图像梯度。SIFT描述子在每个像素点附近选取子区域网格，在每个子区域内计算图像梯度方向直方图。每个子区域的直方图拼接起来组成描述子向量。SIFT描述子的标准设置使用$4\times 4$的子区域，每个子区域使用8个小区间的方向直方图，会产生共128个小区间的直方图$(4\times 4\times 8=128)$。

2.2.3 检测兴趣点

可以使用开源工具包VLFeat提供的二进制文件来计算图像的SIFT特征。用完整的Python实现SIFT特征的所有步骤可能效率不是很高。VLFeat库使用C语言来写的，但是我们可以使用该库提供的命令行接口。由于Python包装器对平台的依赖性，安装Python包装器在某些平台上需要一定的技巧，可以使用二进制文件版本。

2.2.4 匹配描述子

对于讲一幅图像中的特征匹配到另一幅图像的特征，一种稳健的准则（Lowe提出）是使用这两个特征距离和两个最匹配特征距离的比率。相比于图像中的其他特征，该准则保证能够找到足够相似的唯一特征。使用该方法可以使错误的匹配数降低。

2.2.5 代码实现

import cv2
import numpy as np

def cv_show(name, image):
    cv2.imshow(name, image)
    cv2.waitKey(0)
    cv2.destroyAllWindows()

def detectAndDescribe(image):
    # 将彩色图片转换成灰度图
    gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)
    # 建立SIFT生成器
    descriptor = cv2.xfeatures2d.SIFT_create()
    # 检测SIFT特征点，并计算描述子
    (kps, features) = descriptor.detectAndCompute(image, None)
    # 将结果转换成NumPy数组
    kps = np.float32([kp.pt for kp in kps])
    # 返回特征点集，及对应的描述特征
    return (kps, features)

def matchKeypoints(kpsA, kpsB, featuresA, featuresB, ratio = 0.75, reprojThresh = 4.0):
    # 建立暴力匹配器
    matcher = cv2.BFMatcher()
    # 使用KNN检测来自A、B图的SIFT特征匹配对，K=2
    rawMatches = matcher.knnMatch(featuresA, featuresB, 2)
    matches = []
    for m in rawMatches:
        # 当最近距离跟次近距离的比值小于ratio值时，保留此匹配对
        if len(m) == 2 and m[0].distance < m[1].distance * ratio:
        # 存储两个点在featuresA, featuresB中的索引值
            matches.append((m[0].trainIdx, m[0].queryIdx))
        # 当筛选后的匹配对大于4时，计算视角变换矩阵
    if len(matches) > 4:
        # 获取匹配对的点坐标
        ptsA = np.float32([kpsA[i] for (_, i) in matches])
        ptsB = np.float32([kpsB[i] for (i, _) in matches])
        # 计算视角变换矩阵
        (H, status) = cv2.findHomography(ptsA, ptsB, cv2.RANSAC, reprojThresh)
        # 返回结果
        return (matches, H, status)
    # 如果匹配对小于4时，返回None
    return None

def drawMatches(imageA, imageB, kpsA, kpsB, matches, status):
    # 初始化可视化图片，将A、B图左右连接到一起
    (hA, wA) = imageA.shape[:2]
    (hB, wB) = imageB.shape[:2]
    vis = np.zeros((max(hA, hB), wA + wB, 3), dtype="uint8")
    vis[0:hA, 0:wA] = imageA
    vis[0:hB, wA:] = imageB
    # 联合遍历，画出匹配对
    for ((trainIdx, queryIdx), s) in zip(matches, status):
        # 当点对匹配成功时，画到可视化图上
        if s == 1:
            # 画出匹配对
            ptA = (int(kpsA[queryIdx][0]), int(kpsA[queryIdx][1]))
            ptB = (int(kpsB[trainIdx][0]) + wA, int(kpsB[trainIdx][1]))
            cv2.line(vis, ptA, ptB, (0, 255, 0), 1)
    cv_show("drawImg", vis)
    # 返回可视化结果
    return vis

def stitch(imageA,imageB, ratio=0.75, reprojThresh=4.0,showMatches=False):
    #检测A、B图片的SIFT关键特征点，并计算特征描述子
    (kpsA, featuresA) = detectAndDescribe(imageA)
    (kpsB, featuresB) = detectAndDescribe(imageB)
    # 匹配两张图片的所有特征点，返回匹配结果
    M = matchKeypoints(kpsA, kpsB, featuresA, featuresB, ratio, reprojThresh)
    # 如果返回结果为空，没有匹配成功的特征点，退出算法
    if M is None:
        return None
    # 否则，提取匹配结果
    # H是3x3视角变换矩阵
    (matches, H, status) = M
    # 将图片A进行视角变换，result是变换后图片
    result = cv2.warpPerspective(imageA, H, (imageA.shape[1] + imageB.shape[1], imageA.shape[0]))
    cv_show('result', result)
    # 将图片B传入result图片最左端
    result[0:imageB.shape[0], 0:imageB.shape[1]] = imageB
    cv_show('result', result)
    # 检测是否需要显示图片匹配
    if showMatches:
        # 生成匹配图片
        vis = drawMatches(imageA, imageB, kpsA, kpsB, matches, status)
        # 返回结果
        return (result, vis)
    # 返回匹配结果
    return result

# 读取图像
imageA = cv2.imread('/Users/hanjingliu/workspace/daily/course/computer_vision/data/data0/2.pic.jpg')
cv_show("imageA", imageA)
imageB = cv2.imread('/Users/hanjingliu/workspace/daily/course/computer_vision/data/data0/4.pic.jpg')
cv_show("imageB", imageB)
# 计算SIFT特征点和特征向量
(kpsA, featuresA) = detectAndDescribe(imageA)
(kpsB, featuresB) = detectAndDescribe(imageB)
# 基于最近邻和随机取样一致性得到一个单应性矩阵
(matches,H,status) = matchKeypoints(kpsA, kpsB, featuresA, featuresB)
print(H)
# 绘制匹配结果
drawMatches(imageA, imageB, kpsA, kpsB, matches, status)
# 拼接
stitch(imageA, imageB)

代码实现过程：首先读取同一地点的两张不同角度图片，然后进行特征点匹配，最后实现图像拼接。
运行结果如下：
读取到的两张图片：


对读取到的两张图片进行特征匹配：

图像拼接：

2.3 匹配地理标记图像

可以从谷歌提供的照片共享服务Panoramio获得地理标记图像。像许多网络资源一样，Panoramio提供一个API接口，方便用户使用程序访问这些内容。

2.3.1 使用局部描述子匹配

对已下载的图像提取局部描述子。在这种情况下，我们将使用SIFT特征描述子。我们假设已经对这些图像使用SIFT特征提取代码进行了处理，并将特征保存在和图像同名（但文件名后缀是.sift，而不是.jpg）的文件中。

2.3.2 可视化连接的图像

我们首选通过图像间是否具有匹配的局部描述子来定义图像间的连接，然后可视化这些连接情况。为了完成可视化，我们可以在图中显示这些图像，图的边代表连接。我们将会使用pydot工具包，该工具包是功能强大的GraphViz图形库的Python接口。Pydot使用Pyparsing和GraphViz。
为了创建显示可能图像组的图，如果匹配的数目高于一个阈值，我们使用边来连接相应的图像节点。为了得到图中的图像，需要使用图像的全路径。为了使图像看起来更漂亮，我们需要将每幅图像尺度化为缩略图形式。

2.3.3 代码实现

import cv2 as cv
from pylab import *
import os
import pydotplus as pydot

maxsize = (100, 100)  # 定义缩略图的大小
path = r'course/computer_vision/data/all_data'

# 读取整个文件夹的图片
def read_path(pathname):
    imgname_list = os.listdir(pathname)
    img_list = []
    i = 0
    # 图片列表
    for imgname in imgname_list:
        if imgname.endswith('.jpg'):
            img = cv.imread(pathname + '/' + imgname)
            img_n = cv.resize(img, maxsize, cv.INTER_AREA)
            filename = path + str(i) + '.png'
            cv.imwrite(filename, img_n)  # need temporary files of the right size
            i = i + 1
            print(i)
    return img_list

list = read_path(r'course/computer_vision/data/all_data')
print(list)

# 读取整个文件夹的图片
def read_path(pathname):
    imgname_list = os.listdir(pathname)
    img_list = []
    # 图片列表
    for imgname in imgname_list:
        if imgname.endswith('.jpg'):
            img = cv.imread(pathname + '/' + imgname)
            img_list.append(img)
    return img_list

img_list = read_path(r'course/computer_vision/data/all_data')
nbr_images = len(img_list)
match_scores = zeros((nbr_images, nbr_images))

for i in range(nbr_images):
    for j in range(i, nbr_images):  # only compute upper triangle
        print('comparing ', i, j)
        sift = cv.xfeatures2d.SIFT_create()
        kp1, des1 = sift.detectAndCompute(img_list[i], None)
        kp2, des2 = sift.detectAndCompute(img_list[j], None)
        # BFMatch匹配
        bf = cv.BFMatcher(cv.NORM_L2)
        matches = bf.knnMatch(des1, des2, k=2)
        # 储存差距小的优秀匹配点
        goodMatches = []
        for m, n in matches:
            if m.distance < 0.5 * n.distance:
                goodMatches.append(m)
        # 计算优秀匹配点的和
        nbr_matches = len(goodMatches)
        # 向match_scores赋值
        print('number of matches = ', nbr_matches)
        match_scores[i, j] = nbr_matches

# 复制
for i in range(nbr_images):
    for j in range(i + 1, nbr_images):  # no need to copy diagonal 不用复制自我匹配的对角线
        match_scores[j, i] = match_scores[i, j]
# 可视化
threshold = 2  #  至少2个以上匹配点就可以算是联系
g = pydot.Dot(graph_type='graph')  # 不需要有向图
maxsize = (100, 100)  # 定义缩略图的大小
path = r'course/computer_vision/data/all_data'
#两两配对
for i in range(nbr_images):
    for j in range(i + 1, nbr_images):
        if match_scores[i, j] > threshold:
            filename = path + str(i) + '.png'
            g.add_node(pydot.Node(str(i), fontcolor='transparent', shape='rectangle', image=filename))
            filename = path + str(j) + '.png'
            g.add_node(pydot.Node(str(j), fontcolor='transparent', shape='rectangle', image=filename))
            g.add_edge(pydot.Edge(str(i), str(j)))
#绘制S地理标记SIFT匹配图
g.write_jpg('sift.jpg')

运行结果：

==结果分析：==从运行结果中，我们可以看出：三个不同场景拍出的图片却有两个场景中的照片连接在了一起。正常情况下，三个场景的效果应该如上图最右边的一样。这里我们分析得出可能是因为局部描述子匹配方法在两个不同位置的图像中找到了相似的特征点，从而将它们误匹配在一起。想要解决这个问题，可以尝试使用更准确和鲁棒的特征点检测算法来减少错误匹配，并使用更强大的局部描述子来将特征点正确地匹配在一起。另外，也可以考虑使用一些去除错误匹配的方法，例如基于距离、几何一致性和局部一致性的过滤器，以进一步提高匹配的准确性。