import cv2
import numpy as np
import matplotlib
# 强制使用Agg后端，避免GUI问题
matplotlib.use('Agg')
import matplotlib.pyplot as plt
from matplotlib.patches import Polygon
import os
# 设置环境变量，禁用Qt后端
os.environ['QT_QPA_PLATFORM'] = 'offscreen'
os.environ['OPENCV_VIDEOIO_BACKEND'] = 'FFMPEG'
class ImageAlignmentEvaluator:
def __init__(self, img1, img2):
self.img1 = img1
self.img2 = img2
self.H = None
self.F = None
self.homography_score = 0
self.fundamental_score = 0
self.pts1 = None
self.pts2 = None
self.matches = None
def extract_keypoints(self):
"""提取特征点和匹配点"""
# 转换为灰度图像
if len(self.img1.shape) == 3:
gray1 = cv2.cvtColor(self.img1, cv2.COLOR_BGR2GRAY)
gray2 = cv2.cvtColor(self.img2, cv2.COLOR_BGR2GRAY)
else:
gray1 = self.img1
gray2 = self.img2
# 使用ORB特征检测器，避免SIFT的专利问题
orb = cv2.ORB_create(1000)
# 检测关键点和描述符
kp1, des1 = orb.detectAndCompute(gray1, None)
kp2, des2 = orb.detectAndCompute(gray2, None)
if des1 is None or des2 is None:
print("无法提取到足够的特征点，尝试使用SIFT...")
# 回退到SIFT
sift = cv2.SIFT_create()
kp1, des1 = sift.detectAndCompute(gray1, None)
kp2, des2 = sift.detectAndCompute(gray2, None)
if des1 is None or des2 is None:
print("仍然无法提取特征点")
return None, None, None, None, None
# 使用BFMatcher进行匹配
bf = cv2.BFMatcher(cv2.NORM_HAMMING, crossCheck=True)
try:
matches = bf.match(des1, des2)
matches = sorted(matches, key=lambda x: x.distance)
except:
# 如果HAMMING距离不行，尝试L2距离
bf = cv2.BFMatcher()
matches = bf.knnMatch(des1, des2, k=2)
# 应用Lowe's比率测试
good_matches = []
for m, n in matches:
if m.distance < 0.7 * n.distance:
good_matches.append(m)
matches = good_matches
if len(matches) < 10:
print(f"匹配点太少: {len(matches)}")
return None, None, None, None, None
# 获取匹配点的坐标
pts1 = np.float32([kp1[m.queryIdx].pt for m in matches]).reshape(-1, 1, 2)
pts2 = np.float32([kp2[m.trainIdx].pt for m in matches]).reshape(-1, 1, 2)
return pts1, pts2, matches, kp1, kp2
def compute_homography(self, pts1, pts2):
"""计算单应性矩阵H"""
if len(pts1) >= 4:
H, mask = cv2.findHomography(pts2, pts1, cv2.RANSAC, 5.0)
return H, mask
return None, None
def compute_fundamental(self, pts1, pts2):
"""计算基础矩阵F"""
if len(pts1) >= 8:
F, mask = cv2.findFundamentalMat(pts1, pts2, cv2.FM_RANSAC)
return F, mask
return None, None
def evaluate_homography(self, H, pts1, pts2):
"""评估单应性矩阵的质量"""
if H is None:
return 0
try:
# 将pts2变换到图像1的坐标系
pts2_transformed = cv2.perspectiveTransform(pts2, H)
# 计算重投影误差
errors = np.linalg.norm(pts1 - pts2_transformed, axis=2)
mean_error = np.mean(errors)
# 计算内点比例
inlier_ratio = np.sum(errors < 5.0) / len(errors)
# 综合评分（误差越小、内点越多，评分越高）
score = (1 / (1 + mean_error)) * inlier_ratio
return score
except:
return 0
def evaluate_fundamental(self, F, pts1, pts2):
"""评估基础矩阵的质量"""
if F is None:
return 0
try:
# 计算极线约束误差
lines1 = cv2.computeCorrespondEpilines(pts2.reshape(-1, 2), 2, F)
lines1 = lines1.reshape(-1, 3)
errors = []
for (x1, y1), (a, b, c) in zip(pts1.reshape(-1, 2), lines1):
# 点到直线的距离
distance = abs(a*x1 + b*y1 + c) / np.sqrt(a*a + b*b)
errors.append(distance)
mean_error = np.mean(errors)
inlier_ratio = np.sum(np.array(errors) < 1.0) / len(errors)
score = (1 / (1 + mean_error)) * inlier_ratio
return score
except:
return 0
def compare_matrices(self):
"""比较H和F矩阵"""
# 提取特征点
self.pts1, self.pts2, self.matches, kp1, kp2 = self.extract_keypoints()
if self.pts1 is None or len(self.matches) < 10:
print("匹配点太少，无法进行可靠的矩阵计算")
# 尝试使用角点检测作为备选方案
return self.fallback_alignment()
print(f"找到 {len(self.matches)} 个良好的匹配点")
# 计算H矩阵
self.H, H_mask = self.compute_homography(self.pts1, self.pts2)
self.homography_score = self.evaluate_homography(self.H, self.pts1, self.pts2) if self.H is not None else 0
# 计算F矩阵
self.F, F_mask = self.compute_fundamental(self.pts1, self.pts2)
self.fundamental_score = self.evaluate_fundamental(self.F, self.pts1, self.pts2) if self.F is not None else 0
print(f"单应性矩阵(H)评分: {self.homography_score:.4f}")
print(f"基础矩阵(F)评分: {self.fundamental_score:.4f}")
return True
def fallback_alignment(self):
"""备选对齐方案：使用简单的平移对齐"""
print("使用备选对齐方案...")
# 转换为灰度
gray1 = cv2.cvtColor(self.img1, cv2.COLOR_BGR2GRAY)
gray2 = cv2.cvtColor(self.img2, cv2.COLOR_BGR2GRAY)
# 使用模板匹配找到最佳平移
result = cv2.matchTemplate(gray1, gray2, cv2.TM_CCOEFF_NORMED)
_, max_val, _, max_loc = cv2.minMaxLoc(result)
h, w = gray2.shape
x, y = max_loc
# 创建简单的平移矩阵
self.H = np.array([[1, 0, x], [0, 1, y], [0, 0, 1]], dtype=np.float32)
self.homography_score = max_val
self.fundamental_score = 0
print(f"备选方案评分: {self.homography_score:.4f}")
return True
def warp_with_homography(self, alpha=0.5):
"""使用单应性矩阵变换图像"""
if self.H is None:
return None, None, None
h1, w1 = self.img1.shape[:2]
h2, w2 = self.img2.shape[:2]
try:
# 变换图像2到图像1的坐标系
img2_warped = cv2.warpPerspective(self.img2, self.H, (w1, h1))
# 混合图像
blended = cv2.addWeighted(self.img1, 1-alpha, img2_warped, alpha, 0)
# 找到变换后图像的边界框
corners = np.array([[0, 0, 1], [w2, 0, 1], [w2, h2, 1], [0, h2, 1]]).T
transformed_corners = self.H @ corners
transformed_corners = transformed_corners / transformed_corners[2, :]
transformed_corners = transformed_corners[:2, :].T.astype(int)
return blended, img2_warped, transformed_corners
except Exception as e:
print(f"H矩阵变换错误: {e}")
return None, None, None
def warp_with_fundamental(self, alpha=0.5):
"""使用基础矩阵相关的变换"""
if self.F is None or self.pts1 is None:
return None, None, None
try:
# 使用F矩阵筛选内点重新计算H矩阵
lines1 = cv2.computeCorrespondEpilines(self.pts2.reshape(-1, 2), 2, self.F)
lines1 = lines1.reshape(-1, 3)
f_inliers = []
for i, ((x1, y1), (a, b, c)) in enumerate(zip(self.pts1.reshape(-1, 2), lines1)):
distance = abs(a*x1 + b*y1 + c) / np.sqrt(a*a + b*b)
if distance < 1.0: # F矩阵内点阈值
f_inliers.append(i)
if len(f_inliers) >= 4:
# 使用F矩阵的内点重新计算H矩阵
pts1_f = self.pts1[f_inliers]
pts2_f = self.pts2[f_inliers]
H_f, _ = self.compute_homography(pts1_f, pts2_f)
if H_f is not None:
h1, w1 = self.img1.shape[:2]
h2, w2 = self.img2.shape[:2]
img2_warped = cv2.warpPerspective(self.img2, H_f, (w1, h1))
blended = cv2.addWeighted(self.img1, 1-alpha, img2_warped, alpha, 0)
# 计算边界框
corners = np.array([[0, 0, 1], [w2, 0, 1], [w2, h2, 1], [0, h2, 1]]).T
transformed_corners = H_f @ corners
transformed_corners = transformed_corners / transformed_corners[2, :]
transformed_corners = transformed_corners[:2, :].T.astype(int)
return blended, img2_warped, transformed_corners
except Exception as e:
print(f"F矩阵变换错误: {e}")
return None, None, None
def create_comparison_image(self, alpha=0.5):
"""创建H和F矩阵的对比图像"""
blended_h, warped_h, corners_h = self.warp_with_homography(alpha)
blended_f, warped_f, corners_f = self.warp_with_fundamental(alpha)
if blended_h is None:
print("无法使用H矩阵进行变换")
return None, None, None
h, w = self.img1.shape[:2]
# 创建对比图像
comparison_img = np.zeros((h, w*2, 3), dtype=np.uint8)
# 左侧：H矩阵结果
comparison_img[0:h, 0:w] = blended_h
# 右侧：F矩阵结果（如果可用）或原图
if blended_f is not None:
comparison_img[0:h, w:w*2] = blended_f
else:
# 如果F矩阵不可用，显示提示信息
temp_img = self.img1.copy()
cv2.putText(temp_img, "F Matrix Not Available", (w//4, h//2),
cv2.FONT_HERSHEY_SIMPLEX, 1, (0, 0, 255), 2)
comparison_img[0:h, w:w*2] = temp_img
return comparison_img, corners_h, corners_f
def visualize_comprehensive_results(self, alpha=0.5):
"""综合可视化结果"""
try:
# 获取两种方法的变换结果
blended_h, warped_h, corners_h = self.warp_with_homography(alpha)
blended_f, warped_f, corners_f = self.warp_with_fundamental(alpha)
if blended_h is None:
print("无法生成可视化结果")
return
# 创建对比图像
comparison_img, corners_h, corners_f = self.create_comparison_image(alpha)
# 创建综合可视化图表
fig = plt.figure(figsize=(20, 12))
# 1. 原始图像
ax1 = plt.subplot(2, 4, 1)
ax1.imshow(cv2.cvtColor(self.img1, cv2.COLOR_BGR2RGB))
ax1.set_title('Image 1 (Reference)', fontsize=12, fontweight='bold')
ax1.axis('off')
ax2 = plt.subplot(2, 4, 2)
ax2.imshow(cv2.cvtColor(self.img2, cv2.COLOR_BGR2RGB))
ax2.set_title('Image 2 (To Transform)', fontsize=12, fontweight='bold')
ax2.axis('off')
# 2. H矩阵变换结果
ax3 = plt.subplot(2, 4, 3)
ax3.imshow(cv2.cvtColor(blended_h, cv2.COLOR_BGR2RGB))
if corners_h is not None:
poly_h = Polygon(corners_h, fill=False, edgecolor='red', linewidth=3)
ax3.add_patch(poly_h)
ax3.set_title('H Matrix Result\n(Score: {:.3f})'.format(self.homography_score),
fontsize=12, fontweight='bold')
ax3.axis('off')
# 3. F矩阵变换结果
ax4 = plt.subplot(2, 4, 4)
if blended_f is not None:
ax4.imshow(cv2.cvtColor(blended_f, cv2.COLOR_BGR2RGB))
if corners_f is not None:
poly_f = Polygon(corners_f, fill=False, edgecolor='blue', linewidth=3)
ax4.add_patch(poly_f)
ax4.set_title('F Matrix Result\n(Score: {:.3f})'.format(self.fundamental_score),
fontsize=12, fontweight='bold')
else:
ax4.imshow(cv2.cvtColor(self.img1, cv2.COLOR_BGR2RGB))
ax4.set_title('F Matrix Not Available', fontsize=12, fontweight='bold')
ax4.axis('off')
# 4. 对比图像
ax5 = plt.subplot(2, 2, 3)
if comparison_img is not None:
ax5.imshow(cv2.cvtColor(comparison_img, cv2.COLOR_BGR2RGB))
# 添加分割线
h, w = self.img1.shape[:2]
ax5.axvline(x=w, color='white', linewidth=3, linestyle='--')
ax5.text(w//4, 30, 'H Matrix', fontsize=14, fontweight='bold',
color='red', ha='center', backgroundcolor='white')
ax5.text(w*1.5, 30, 'F Matrix', fontsize=14, fontweight='bold',
color='blue', ha='center', backgroundcolor='white')
ax5.set_title('Method Comparison', fontsize=14, fontweight='bold')
ax5.axis('off')
# 5. 评分比较
ax6 = plt.subplot(2, 2, 4)
matrices = ['H Matrix', 'F Matrix']
scores = [self.homography_score, self.fundamental_score]
colors = ['lightcoral', 'lightblue']
bars = ax6.bar(matrices, scores, color=colors, alpha=0.7)
ax6.set_ylabel('Score', fontsize=12)
ax6.set_title('Matrix Score Comparison', fontsize=14, fontweight='bold')
# 添加数值标签
for bar, score in zip(bars, scores):
height = bar.get_height()
ax6.text(bar.get_x() + bar.get_width()/2., height + 0.01,
f'{score:.3f}', ha='center', va='bottom', fontsize=12)
plt.tight_layout()
plt.savefig('comprehensive_comparison_result.png', dpi=300, bbox_inches='tight')
print("可视化结果已保存为 'comprehensive_comparison_result.png'")
# 保存对比图像
if comparison_img is not None:
cv2.imwrite('method_comparison.jpg', comparison_img)
cv2.imwrite('H_matrix_result.jpg', blended_h)
if blended_f is not None:
cv2.imwrite('F_matrix_result.jpg', blended_f)
print("结果已保存:")
print("- comprehensive_comparison_result.png")
print("- method_comparison.jpg")
print("- H_matrix_result.jpg")
if blended_f is not None:
print("- F_matrix_result.jpg")
# 不要显示图表，直接保存
plt.close(fig)
except Exception as e:
print(f"可视化错误: {e}")
import traceback
traceback.print_exc()
import os
def main():
# 设置环境变量
os.environ['QT_QPA_PLATFORM'] = 'offscreen'
# 读取图像
path_="/home/dongdong/2project/0data/RTK/data_4_city/300_locatiopn_2pm/images"
# 读取图像
img1_path = os.path.join(path_, 'DJI_0528.JPG')
img2_path = os.path.join(path_, 'DJI_0559.JPG')
try:
print("正在读取图像...")
img1 = cv2.imread(img1_path)
img2 = cv2.imread(img2_path)
if img1 is None or img2 is None:
print("错误: 无法读取图像，请检查文件路径")
print(f"图像1路径: {img1_path}")
print(f"图像2路径: {img2_path}")
return
print(f"图像1尺寸: {img1.shape}")
print(f"图像2尺寸: {img2.shape}")
# 调整图像大小以保持一致性
h1, w1 = img1.shape[:2]
h2, w2 = img2.shape[:2]
if h1 != h2 or w1 != w2:
print("图像尺寸不一致，进行调整...")
target_size = (min(w1, w2), min(h1, h2))
img1 = cv2.resize(img1, target_size)
img2 = cv2.resize(img2, target_size)
print(f"调整后尺寸: {img1.shape}")
# 创建评估器实例
evaluator = ImageAlignmentEvaluator(img1, img2)
# 比较矩阵
print("正在计算特征点和矩阵...")
if not evaluator.compare_matrices():
print("矩阵计算失败，使用备选方案")
# 生成综合可视化结果
print("正在生成可视化结果...")
evaluator.visualize_comprehensive_results(alpha=0.5)
print("处理完成！")
except Exception as e:
print(f"发生错误: {e}")
import traceback
traceback.print_exc()
if __name__ == "__main__":
main()
浙公网安备 33010602011771号