运用视觉里程计尝试观测月球：ORB+对极几何+三角化+尺度重建

一、背景

前两年听了一场关于单目多视角Slam技术的讲座，尽管精度稳定性都不如双目或雷达，计算起来也更加复杂，但是硬件简化和研究论文效果演示还是给我留下了深刻印象。

也买了《视觉SLAM十四讲》，高博士写的很好，里面的例子循序渐进，这里也多有参考引用里面的内容。

最近在用天文望远镜看月亮的时候我就在想，当我连续观测月亮，地球在自转，根据SFM（运动恢复结构）思想，用两个位置（由地球自转产生）模拟一个双目系统，本质上应该是是单目多视角系统，但是忽略了时间序列上的信息积累和全局优化，比如BA优化（Bundle Adjustment，光束法平差），简化成“特征匹配+几何求解”。

这个“双目”的基线相对于地月系统太短，视差太小，只是做着玩玩试试看。

先上素材，来自大气的扰动还是很剧烈的。

二、结论

先说结论：这是一次失败！！！的尝试，只是一次有趣的学习，为一时天马行空的想法证明了错误路线。尽管SFM在天文观测中有运用，但是主要是环绕观测。

三、关键知识点

1. ORB特征检测与匹配技术

1.1 ORB算法原理

ORB (Oriented FAST and Rotated BRIEF) 是一种高效的特征检测与描述算法，特别适合在资源受限的环境中使用。

FAST角点检测

• 检测图像中亮度变化显著的像素点
• 通过圆形邻域比较快速识别角点
• 计算复杂度低，满足实时性要求

BRIEF描述子

• 生成二进制字符串描述特征点周围区域
• 比传统浮点描述子(SIFT/SURF)计算更快
• 存储空间需求小，匹配效率高

1.2 环境适应性改进

• 自适应阈值调整应对极端光照变化
• 多尺度特征检测适应不同分辨率需求
• 抗噪声设计提高在尘埃环境下的鲁棒性

2. 对极几何理论与应用

2.1 几何基础

对极几何描述了两个相机视图间的几何关系，是立体视觉的基础理论。

2.2 关键概念

• 基线: 两个相机光心的连线
• 对极平面: 包含基线和空间点的平面
• 对极点: 相机光心在对方图像平面上的投影

2.3 数学描述

基础矩阵(Fundamental Matrix)

F = K₂⁻ᵀEK₁⁻¹

其中E为本质矩阵，K₁,K₂为相机内参矩阵。

约束方程

x₂ᵀFx₁ = 0

该方程表达了对应点必须满足的几何约束。

3. 三角化重建原理

3.1 DLT算法

直接线性变换(Direct Linear Transform)是一种常用的三角化方法。

3.2 数学推导

给定两个视图的投影矩阵P₁,P₂和对应点x₁,x₂，空间点X可通过求解齐次线性方程组得到。

4. 尺度重建与优化

4.1 尺度不确定性问题

单目视觉系统存在固有的尺度模糊性，需要额外信息确定真实尺度。

4.2 解决方案

IMU融合

• 利用惯性测量单元提供运动信息
• 通过传感器融合提高位姿估计精度
• 补偿视觉观测的短期漂移

已知参考物

• 利用已知尺寸的地标进行尺度校正（这里简单演示了利用已知月球直径）
• 结合轨迹信息优化尺度因子
• 多视图一致性约束提高重建可靠性

5. 相机内参

CMOS型号：SONY IMX335 1/2.8
1/2.8 英寸传感器：
– 对角线 = 16 / 2.8 ≈ 5.71 mm
– 宽 = 5.704 mm
– 高 = 4.28 mm
– 宽高比 = 4:3
分辨率：2592*1944
镜头焦距：500mm
超长焦相机，难以标定，估计的参数
┌ ┐
│ fx 0 cx │
K = │ 0 fy cy │
│ 0 0 1 │
└ ┘
fx = f / dx = f × width / sensor_width
fy = f / dy =f × height / sensor_height
dx和 dy： 是相机传感器（如CMOS）上单个像素的物理尺寸（单位：毫米/像素）。
| `K(0,0)` | | `fx` | **x 方向焦距** | 像素 |
| `K(0,1)` | 0 | – | 倾斜系数（通常为 0） | – |
| `K(0,2)` | 1296 | `cx` | **主点 x 坐标**（光心水平位置） | 像素 |
| `K(1,0)` | 0 | – | 倾斜系数（通常为 0） | – |
| `K(1,1)` | | `fy` | **y 方向焦距** | 像素 |
| `K(1,2)` | 972 | `cy` | **主点 y 坐标**（光心垂直位置） | 像素 |
| `K(2,0)` | 0 | – | 齐次坐标固定值 | – |
| `K(2,1)` | 0 | – | 齐次坐标固定值 | – |
| `K(2,2)` | 1 | – | 齐次坐标固定值 | – |

四、软件流程图

五、C++实现

ORB_pattern 是 Willow Garage 团队在 2011 年 通过机器学习从 PASCAL VOC 等数据集 训练得到的标准 ORB 描述子采样模式，属于计算机视觉领域的公开算法参数，已集成到 OpenCV 等主流库中，无需自行训练。

// 文件: visual_odom.cpp
// 描述: 基于单目视觉的月球视频处理程序，用于特征匹配、运动估计、三角测量及尺度恢复。
// 核心算法: ORB特征提取与匹配，对极几何，RANSAC，三角测量。

#include <opencv2/opencv.hpp>
#include <chrono>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <limits>
#include <string>
#include <sys/stat.h>
#include <unistd.h>
#include <vector>

// SSE intrinsics for efficient Hamming distance calculation
#ifdef __SSE4_2__
#include <nmmintrin.h>
#include <popcntintrin.h>
#endif

// 使用描述符类型表示一个256位的ORB描述符，由8个32位无符号整数组成。
typedef std::vector<uint32_t> OrbDescriptor;

// 全局常量定义
namespace constants
{
// 月球物理参数 (单位: km)
const double kMoonDiameterPixel = 2352.17; // 月球在图像中的像素直径
const double kScaleKmPerPixel = 1.478; // 比例尺 (km/pixel)，由标定得到
const double kMoonDistanceKm = 384400.0; // 平均地月距离
const double kMoonRadiusKm = 1738.14; // 月球平均半径
const double kMoonDiameterKm = 3476.28; // 月球平均直径
const double kEarthRotationSpeed = 0.403; // 地球自转线速度 (km/s) @ 北纬30度
} // namespace constants

// 外部定义的ORB点对比较模式表，用于计算旋转不变的BRIEF描述子。
// 格式: 对于256位描述子，共有256个比较对，每个点对由4个int定义(p.x, p.y, q.x, q.y)。
extern int ORB_pattern[256 * 4];

// 函数声明
void ComputeOrbDescriptors(const cv::Mat &grayscale_image,
std::vector<cv::KeyPoint> *keypoints,
std::vector<OrbDescriptor> *descriptors);

void BruteForceMatch(const std::vector<OrbDescriptor> &descriptors1,
const std::vector<OrbDescriptor> &descriptors2,
std::vector<cv::DMatch> *matches);

cv::Point2d PixelToNormalizedPlane(const cv::Point2d &pixel_point,
const cv::Mat &camera_matrix);

void EstimatePoseFrom2d2dMatches(
const std::vector<cv::KeyPoint> &keypoints1,
const std::vector<cv::KeyPoint> &keypoints2,
const std::vector<cv::DMatch> &matches,
const cv::Mat &camera_matrix,
cv::Mat *rotation_matrix,
cv::Mat *translation_vector);

void TriangulatePoints(
const std::vector<cv::KeyPoint> &keypoints1,
const std::vector<cv::KeyPoint> &keypoints2,
const std::vector<cv::DMatch> &matches,
const cv::Mat &camera_matrix,
const cv::Mat &rotation_matrix,
const cv::Mat &translation_vector,
std::vector<cv::Point3d> *triangulated_points);

inline cv::Scalar GetDepthColor(float depth,
float upper_threshold = 50.0f,
float lower_threshold = 10.0f);

void FilterMatchesByDepthConsistency(
const std::vector<cv::Point3d> &points,
const std::vector<cv::DMatch> &input_matches,
std::vector<cv::DMatch> *filtered_matches,
double depth_threshold_factor = 3.0);

cv::Point3d FindFarthestPointPairRansac(
const std::vector<cv::Point3d> &points,
int *best_index1,
int *best_index2,
int iterations = 1000,
double inlier_threshold = 0.1);

bool ValidateSphereAssumption(const std::vector<cv::Point3d> &points_km);

// 主函数
int main(int argc, char **argv)
{
// 视频文件路径 (可考虑通过命令行参数传入)
const std::string kVideoFilePath = "../moon_video.mp4";
const std::string kOutputDir = "matches_frame";

// 创建输出目录
if (access(kOutputDir.c_str(), F_OK) != 0)
{
if (mkdir(kOutputDir.c_str(), 0755) == 0)
{
std::cout << "创建输出文件夹: " << kOutputDir << std::endl;
}
else
{
std::cerr << "错误: 无法创建输出文件夹 " << kOutputDir << std::endl;
return -1;
}
}

// 打开视频文件
cv::VideoCapture video_cap(kVideoFilePath);
if (!video_cap.isOpened())
{
std::cerr << "错误: 无法打开视频文件: " << kVideoFilePath << std::endl;
return -1;
}

// 获取视频基本信息
const int kTotalFrameCount =
static_cast<int>(video_cap.get(cv::CAP_PROP_FRAME_COUNT));
const double kFps = video_cap.get(cv::CAP_PROP_FPS);
const int kFrameWidth =
static_cast<int>(video_cap.get(cv::CAP_PROP_FRAME_WIDTH));
const int kFrameHeight =
static_cast<int>(video_cap.get(cv::CAP_PROP_FRAME_HEIGHT));

// 输出视频信息
std::cout << "\\n========== 视频文件信息 ==========" << std::endl;
std::cout << "文件路径: " << kVideoFilePath << std::endl;
std::cout << "总帧数: " << kTotalFrameCount << " 帧" << std::endl;
std::cout << "帧率: " << kFps << " FPS" << std::endl;
std::cout << "分辨率: " << kFrameWidth << " x " << kFrameHeight << " 像素"
<< std::endl;
std::cout << "预计播放时长: "
<< (kTotalFrameCount > 0 && kFps > 0
? std::to_string(kTotalFrameCount / kFps) + " 秒"
: "未知")
<< std::endl;
std::cout << "========================================\\n"
<< std::endl;

// 相机内参矩阵估计
// 传感器: SONY IMX335 1/2.8英寸, 镜头焦距: 500mm
const double kFocalLengthMm = 500.0; // 焦距 (mm)
const double kSensorWidthMm = 5.704; // 传感器宽度 (mm)
const double kSensorHeightMm = 4.28; // 传感器高度 (mm)
const double kFx = kFocalLengthMm * kFrameWidth / kSensorWidthMm;
const double kFy = kFocalLengthMm * kFrameHeight / kSensorHeightMm;
const double kCx = kFrameWidth / 2.0;
const double kCy = kFrameHeight / 2.0;

cv::Mat camera_matrix = (cv::Mat_<double>(3, 3) << kFx, 0, kCx, 0, kFy, kCy,
0, 0, 1);
std::cout << "相机内参矩阵 K:\\n"
<< camera_matrix << std::endl;

// 重新打开两个视频流，用于处理前后帧对
video_cap.release();
cv::VideoCapture cap_front(kVideoFilePath);
cv::VideoCapture cap_back(kVideoFilePath);
if (!cap_front.isOpened() || !cap_back.isOpened())
{
std::cerr << "错误: 重新打开视频文件失败" << std::endl;
return -1;
}

// 处理配置: 帧间隔
const int kFrameInterval = 2;
const double kTimeInterval = kFrameInterval / kFps;
if (!cap_back.set(cv::CAP_PROP_POS_FRAMES, kFrameInterval))
{
std::cerr << "警告: 无法将后半部分视频跳转到第 " << kFrameInterval
<< " 帧" << std::endl;
}

std::cout << "\\n========== 视频处理配置 ==========" << std::endl;
std::cout << "帧间隔: " << kFrameInterval << " 帧" << std::endl;
std::cout << "时间间隔: " << kTimeInterval << " 秒" << std::endl;
std::cout << "预计处理帧对: " << kTotalFrameCount – kFrameInterval << " 对"
<< std::endl;
std::cout << "========================================\\n"
<< std::endl;

// 状态变量初始化
int current_frame_index = 0;
int total_match_count = 0;
double min_estimated_distance = std::numeric_limits<double>::max();
double max_estimated_distance = std::numeric_limits<double>::lowest();
double best_estimated_distance = 0.0;
double best_distance_error_ratio = 1.0;
double best_camera_baseline_km = 0.0;
double best_camera_speed_kms = 0.0;

// 性能计时
auto fps_start_time = std::chrono::steady_clock::now();
int fps_frame_counter = 0;
double current_fps = 0.0;

// 主处理循环
while (true)
{
cv::Mat front_frame_color, back_frame_color;
cv::Mat front_frame_gray, back_frame_gray;

// 读取前后帧
if (!cap_front.read(front_frame_color) ||
!cap_back.read(back_frame_color))
{
std::cout << "视频读取完毕或到达结尾。" << std::endl;
break;
}

// 转换为灰度图
cv::cvtColor(front_frame_color, front_frame_gray, cv::COLOR_BGR2GRAY);
cv::cvtColor(back_frame_color, back_frame_gray, cv::COLOR_BGR2GRAY);

std::cout << "\\n===== 处理帧对 " << (current_frame_index + 1) << " ====="
<< std::endl;
std::cout << "前帧索引: " << current_frame_index
<< ", 后帧索引: " << (current_frame_index + kFrameInterval)
<< std::endl;

// ==================== 第一步: 特征提取 ====================
std::cout << "\\n[步骤1] ORB特征提取…" << std::endl;
auto time_start = std::chrono::steady_clock::now();

std::vector<cv::KeyPoint> keypoints_front, keypoints_back;
std::vector<OrbDescriptor> descriptors_front, descriptors_back;
cv::FAST(front_frame_gray, keypoints_front, 40);
cv::FAST(back_frame_gray, keypoints_back, 40);
ComputeOrbDescriptors(front_frame_gray, &keypoints_front,
&descriptors_front);
ComputeOrbDescriptors(back_frame_gray, &keypoints_back,
&descriptors_back);

auto time_end = std::chrono::steady_clock::now();
std::chrono::duration<double> elapsed = time_end – time_start;
std::cout << "耗时: " << elapsed.count() << " 秒" << std::endl;

// ==================== 第二步: 特征匹配 ====================
std::cout << "\\n[步骤2] 特征匹配 (暴力匹配)…" << std::endl;
time_start = std::chrono::steady_clock::now();

std::vector<cv::DMatch> initial_matches;
BruteForceMatch(descriptors_front, descriptors_back, &initial_matches);
total_match_count += initial_matches.size();

time_end = std::chrono::steady_clock::now();
elapsed = time_end – time_start;
std::cout << "耗时: " << elapsed.count() << " 秒" << std::endl;
std::cout << "本帧匹配对数: " << initial_matches.size() << std::endl;

// ==================== 第三步: 运动估计 (对极几何) ====================
std::cout << "\\n[步骤3] 基于对极几何估计相机运动…" << std::endl;
bool is_pose_estimated = false;
cv::Mat rotation_matrix, translation_vector;
std::vector<cv::DMatch> filtered_matches_by_pose;

if (initial_matches.size() >= 8)
{
EstimatePoseFrom2d2dMatches(keypoints_front, keypoints_back,
initial_matches, camera_matrix,
&rotation_matrix, &translation_vector);
is_pose_estimated = true;
std::cout << "旋转矩阵 R:\\n"
<< rotation_matrix << std::endl;
std::cout << "平移向量 t (单位向量):\\n"
<< translation_vector
<< std::endl;
}
else
{
std::cout << "匹配点对不足 (" << initial_matches.size()
<< ")，跳过运动估计。" << std::endl;
}

// ==================== 第四步: 三角测量 ====================
cv::Mat depth_viz_front = front_frame_color.clone();
cv::Mat depth_viz_back = back_frame_color.clone();
std::vector<cv::Point3d> triangulated_points_normalized;
double current_physical_distance_km = 0.0;
std::vector<cv::DMatch> filtered_matches_by_depth;

if (is_pose_estimated)
{
std::cout << "\\n[步骤4] 三角测量恢复3D点…" << std::endl;
TriangulatePoints(keypoints_front, keypoints_back, initial_matches,
camera_matrix, rotation_matrix,
translation_vector,
&triangulated_points_normalized);

// 基于深度一致性筛选异常匹配点
FilterMatchesByDepthConsistency(triangulated_points_normalized,
initial_matches,
&filtered_matches_by_depth, 3.0);
std::cout << "匹配点数量: " << initial_matches.size() << " -> "
<< filtered_matches_by_depth.size() << " (筛选后)"
<< std::endl;

if (!triangulated_points_normalized.empty())
{
// 计算深度统计并可视化
double min_depth_norm = std::numeric_limits<double>::max();
double max_depth_norm = std::numeric_limits<double>::lowest();
for (const auto &p : triangulated_points_normalized)
{
if (p.z < min_depth_norm)
min_depth_norm = p.z;
if (p.z > max_depth_norm)
max_depth_norm = p.z;
}

for (size_t i = 0; i < filtered_matches_by_depth.size(); ++i)
{
const cv::DMatch &match = filtered_matches_by_depth[i];
float depth1 = triangulated_points_normalized[i].z;
cv::circle(depth_viz_front,
keypoints_front[match.queryIdx].pt, 4,
GetDepthColor(depth1, max_depth_norm,
min_depth_norm),
4);
cv::Mat pt_in_cam2 =
rotation_matrix *
(cv::Mat_<double>(3, 1)
<< triangulated_points_normalized[i].x,
triangulated_points_normalized[i].y,
triangulated_points_normalized[i].z) +
translation_vector;
float depth2 = pt_in_cam2.at<double>(2, 0);
cv::circle(depth_viz_back, keypoints_back[match.trainIdx].pt,
4,
GetDepthColor(depth2, max_depth_norm,
min_depth_norm),
4);
}
}
}

// ==================== 第五步: 尺度恢复 ====================
if (triangulated_points_normalized.size() >= 2)
{
std::cout << "\\n[步骤5] 尺度恢复与物理距离估计…" << std::endl;

// 1. 使用RANSAC寻找最远点对 (假设为月球直径)
int idx1 = 0, idx2 = 0;
FindFarthestPointPairRansac(triangulated_points_normalized, &idx1,
&idx2, 1000, 0.1);
double diameter_norm = cv::norm(
triangulated_points_normalized[idx1] –
triangulated_points_normalized[idx2]);

// 2. 计算尺度因子 (基于已知月球直径约束)
double moon_diameter_norm =
constants::kMoonDiameterPixel / kFx;
double scale_factor =
constants::kMoonDiameterKm / moon_diameter_norm;
std::cout << "计算得到的尺度因子: " << scale_factor << " km/unit"
<< std::endl;

// 3. 恢复物理尺度的平移向量和3D点云
cv::Mat translation_km = translation_vector * scale_factor;
double physical_baseline_km = cv::norm(translation_km);
double physical_camera_speed_kms =
physical_baseline_km / kTimeInterval;
std::cout << "相机间物理基线: " << physical_baseline_km << " km"
<< std::endl;
std::cout << "推算相机运动速度: " << physical_camera_speed_kms
<< " km/s" << std::endl;

std::vector<cv::Point3d> points_km;
for (const auto &p : triangulated_points_normalized)
{
points_km.push_back(p * scale_factor);
}

// 4. 估计相机到月球表面的距离 (最近点深度 + 月球半径)
double min_depth_km = std::numeric_limits<double>::max();
for (const auto &p : points_km)
{
if (p.z < min_depth_km)
min_depth_km = p.z;
}
current_physical_distance_km =
min_depth_km + constants::kMoonRadiusKm;

// 5. 更新最佳/最差估计
if (current_physical_distance_km > max_estimated_distance)
max_estimated_distance = current_physical_distance_km;
if (current_physical_distance_km < min_estimated_distance)
min_estimated_distance = current_physical_distance_km;

double current_error_ratio =
std::fabs(current_physical_distance_km –
constants::kMoonDistanceKm) /
constants::kMoonDistanceKm;
if (current_error_ratio < best_distance_error_ratio)
{
best_distance_error_ratio = current_error_ratio;
best_estimated_distance = current_physical_distance_km;
best_camera_baseline_km = physical_baseline_km;
best_camera_speed_kms = physical_camera_speed_kms;
}

std::cout << "本帧估计地月距离: " << current_physical_distance_km
<< " km" << std::endl;
std::cout << "当前最佳估计: " << best_estimated_distance << " km (误差: "
<< best_distance_error_ratio * 100 << "%)" << std::endl;

// 6. 验证点云是否近似分布在球面上
bool is_sphere_valid = ValidateSphereAssumption(points_km);
if (is_sphere_valid)
{
std::cout << "验证通过: 点云近似分布在球面上。" << std::endl;
}
else
{
std::cout << "警告: 点云分布不均匀，直径假设可能不成立。"
<< std::endl;
}
}

// ==================== 第六步: 可视化与输出 ====================
std::cout << "\\n[步骤6] 生成可视化结果…" << std::endl;
cv::Mat matches_image;
std::vector<cv::DMatch> matches_to_draw =
is_pose_estimated ? filtered_matches_by_depth : initial_matches;
cv::drawMatches(front_frame_color, keypoints_front, back_frame_color,
keypoints_back, matches_to_draw, matches_image);

cv::Mat depth_image;
cv::hconcat(depth_viz_front, depth_viz_back, depth_image);

cv::Mat final_display;
cv::vconcat(matches_image, depth_image, final_display);
cv::line(final_display, cv::Point(0, matches_image.rows),
cv::Point(final_display.cols, matches_image.rows),
cv::Scalar(255, 255, 255), 2);

// 添加状态信息文本
const int kFontFace = cv::FONT_HERSHEY_SIMPLEX;
const double kFontScale = 2.5;
const int kFontThickness = 4;
int line_height = 100;
int y_pos = 100;

auto PutText = [&](const std::string &text)
{
cv::putText(final_display, text, cv::Point(10, y_pos), kFontFace,
kFontScale, cv::Scalar(0, 255, 0), kFontThickness);
y_pos += line_height;
};

PutText("Frame: " + std::to_string(current_frame_index) + " | Interval: " +
std::to_string(kTimeInterval) + " s");
PutText("Matches: " + std::to_string(matches_to_draw.size()));
PutText("Curr Dist: " + std::to_string(current_physical_distance_km) +
" km");
PutText("Best Dist: " + std::to_string(best_estimated_distance) + " km (" +
std::to_string(best_distance_error_ratio * 100) + "%)");
PutText("Cam Speed: " + std::to_string(best_camera_speed_kms) + " km/s");

// 计算并显示实时处理FPS
fps_frame_counter++;
auto now = std::chrono::steady_clock::now();
std::chrono::duration<double> fps_elapsed = now – fps_start_time;
if (fps_frame_counter >= 30)
{
current_fps = fps_frame_counter / fps_elapsed.count();
fps_start_time = now;
fps_frame_counter = 0;
}
std::string fps_text = "Proc FPS: " + std::to_string(static_cast<int>(current_fps));
cv::putText(final_display, fps_text,
cv::Point(final_display.cols – 600, 100), kFontFace,
kFontScale, cv::Scalar(0, 255, 255), kFontThickness);

// 显示结果 (适应屏幕缩放)
const int kMaxDisplayWidth = 1920;
const int kMaxDisplayHeight = 1080;
if (final_display.cols > kMaxDisplayWidth ||
final_display.rows > kMaxDisplayHeight)
{
double scale = std::min(
static_cast<double>(kMaxDisplayWidth) / final_display.cols,
static_cast<double>(kMaxDisplayHeight) / final_display.rows);
cv::Mat resized_display;
cv::resize(final_display, resized_display, cv::Size(), scale, scale,
cv::INTER_LINEAR);
cv::imshow("Visual Odometry for Moon Distance Estimation",
resized_display);
}
else
{
cv::imshow("Visual Odometry for Moon Distance Estimation",
final_display);
}

// 保存结果 (可选)
// std::string save_path = kOutputDir + "/frame_" +
// std::to_string(current_frame_index) + ".png";
// cv::imwrite(save_path, final_display);

// 处理用户按键
char key = cv::waitKey(1);
if (key == 'q' || key == 27)
{ // 'q' 或 ESC
std::cout << "用户中断处理。" << std::endl;
break;
}

current_frame_index++;
}

// ==================== 处理结果摘要 ====================
std::cout << "\\n========== 处理完成 ==========" << std::endl;
std::cout << "处理总帧对: " << current_frame_index << std::endl;
std::cout << "总匹配对数: " << total_match_count << std::endl;
if (current_frame_index > 0)
{
std::cout << "平均每帧匹配对数: "
<< total_match_count / current_frame_index << std::endl;
}
std::cout << "地月距离估计范围: " << min_estimated_distance << " km 到 "
<< max_estimated_distance << " km" << std::endl;
std::cout << "最佳地月距离估计: " << best_estimated_distance << " km (误差: "
<< best_distance_error_ratio * 100 << "%)" << std::endl;
std::cout << "对应相机运动速度: " << best_camera_speed_kms << " km/s"
<< std::endl;
std::cout << "==============================\\n"
<< std::endl;

std::cout << "按任意键关闭窗口…" << std::endl;
cv::waitKey(0);

cv::destroyAllWindows();
return 0;
}

// ============================================================================
// 函数定义
// ============================================================================

void ComputeOrbDescriptors(const cv::Mat &grayscale_image,
std::vector<cv::KeyPoint> *keypoints,
std::vector<OrbDescriptor> *descriptors)
{
const int kHalfPatchSize = 8;
const int kHalfBoundary = 16;
int bad_points = 0;
descriptors->clear();

for (auto &kp : *keypoints)
{
// 检查关键点是否在有效图像区域内
if (kp.pt.x < kHalfBoundary || kp.pt.y < kHalfBoundary ||
kp.pt.x >= grayscale_image.cols – kHalfBoundary ||
kp.pt.y >= grayscale_image.rows – kHalfBoundary)
{
bad_points++;
descriptors->push_back(OrbDescriptor()); // 空描述符
continue;
}

// 计算关键点方向 (图像矩)
float m01 = 0.0f, m10 = 0.0f;
for (int dx = -kHalfPatchSize; dx < kHalfPatchSize; ++dx)
{
for (int dy = -kHalfPatchSize; dy < kHalfPatchSize; ++dy)
{
uchar pixel =
grayscale_image.at<uchar>(kp.pt.y + dy, kp.pt.x + dx);
m10 += dx * pixel;
m01 += dy * pixel;
}
}
float norm = sqrt(m01 * m01 + m10 * m10) + 1e-18f;
float sin_theta = m01 / norm;
float cos_theta = m10 / norm;

// 计算旋转不变的BRIEF描述子 (ORB)
OrbDescriptor desc(8, 0);
for (int i = 0; i < 8; ++i)
{
uint32_t d = 0;
for (int k = 0; k < 32; ++k)
{
int idx = i * 32 + k;
cv::Point2f p(ORB_pattern[idx * 4], ORB_pattern[idx * 4 + 1]);
cv::Point2f q(ORB_pattern[idx * 4 + 2], ORB_pattern[idx * 4 + 3]);

// 旋转点坐标
cv::Point2f rotated_p(cos_theta * p.x – sin_theta * p.y + kp.pt.x,
sin_theta * p.x + cos_theta * p.y + kp.pt.y);
cv::Point2f rotated_q(cos_theta * q.x – sin_theta * q.y + kp.pt.x,
sin_theta * q.x + cos_theta * q.y + kp.pt.y);

// 比较像素强度并设置位
if (grayscale_image.at<uchar>(rotated_p.y, rotated_p.x) <
grayscale_image.at<uchar>(rotated_q.y, rotated_q.x))
{
d |= (1 << k);
}
}
desc[i] = d;
}
descriptors->push_back(desc);
}
std::cout << "忽略的关键点 (靠近边界): " << bad_points << "/"
<< keypoints->size() << std::endl;
}

void BruteForceMatch(const std::vector<OrbDescriptor> &descriptors1,
const std::vector<OrbDescriptor> &descriptors2,
std::vector<cv::DMatch> *matches)
{
const int kMaxHammingDistance = 40;
matches->clear();

for (size_t i1 = 0; i1 < descriptors1.size(); ++i1)
{
if (descriptors1[i1].empty())
continue;

cv::DMatch best_match;
best_match.queryIdx = static_cast<int>(i1);
best_match.trainIdx = -1;
best_match.distance = 256; // 初始化为最大可能距离

for (size_t i2 = 0; i2 < descriptors2.size(); ++i2)
{
if (descriptors2[i2].empty())
continue;

int hamming_distance = 0;
for (int k = 0; k < 8; ++k)
{
#ifdef __SSE4_2__
hamming_distance +=
_mm_popcnt_u32(descriptors1[i1][k] ^ descriptors2[i2][k]);
#else
// 兼容性回退方案
uint32_t xor_val = descriptors1[i1][k] ^ descriptors2[i2][k];
hamming_distance += __builtin_popcount(xor_val);
#endif
}

if (hamming_distance < kMaxHammingDistance &&
hamming_distance < best_match.distance)
{
best_match.distance = hamming_distance;
best_match.trainIdx = static_cast<int>(i2);
}
}

if (best_match.trainIdx != -1)
{
matches->push_back(best_match);
}
}
}

cv::Point2d PixelToNormalizedPlane(const cv::Point2d &pixel_point,
const cv::Mat &camera_matrix)
{
double x = (pixel_point.x – camera_matrix.at<double>(0, 2)) /
camera_matrix.at<double>(0, 0);
double y = (pixel_point.y – camera_matrix.at<double>(1, 2)) /
camera_matrix.at<double>(1, 1);
return cv::Point2d(x, y);
}

void EstimatePoseFrom2d2dMatches(
const std::vector<cv::KeyPoint> &keypoints1,
const std::vector<cv::KeyPoint> &keypoints2,
const std::vector<cv::DMatch> &matches,
const cv::Mat &camera_matrix,
cv::Mat *rotation_matrix,
cv::Mat *translation_vector)
{
// 提取匹配点对坐标
std::vector<cv::Point2f> points1, points2;
for (const auto &match : matches)
{
points1.push_back(keypoints1[match.queryIdx].pt);
points2.push_back(keypoints2[match.trainIdx].pt);
}

// 计算本质矩阵 (使用RANSAC去除误匹配)
cv::Point2d principal_point(camera_matrix.at<double>(0, 2),
camera_matrix.at<double>(1, 2));
double focal_length = camera_matrix.at<double>(0, 0);
cv::Mat essential_matrix;
cv::Mat inliers_mask;
essential_matrix = cv::findEssentialMat(
points1, points2, focal_length, principal_point, cv::RANSAC, 0.999, 1.0,
inliers_mask);

// 从本质矩阵恢复旋转和平移
std::vector<cv::Point2f> points1_inlier, points2_inlier;
for (int i = 0; i < inliers_mask.rows; ++i)
{
if (inliers_mask.at<uchar>(i))
{
points1_inlier.push_back(points1[i]);
points2_inlier.push_back(points2[i]);
}
}

if (points1_inlier.size() >= 5)
{
cv::recoverPose(essential_matrix, points1_inlier, points2_inlier,
*rotation_matrix, *translation_vector, focal_length,
principal_point);
}
else
{
std::cout << "警告: 内点数量不足，姿态估计可能不可靠。" << std::endl;
}
}

void TriangulatePoints(
const std::vector<cv::KeyPoint> &keypoints1,
const std::vector<cv::KeyPoint> &keypoints2,
const std::vector<cv::DMatch> &matches,
const cv::Mat &camera_matrix,
const cv::Mat &rotation_matrix,
const cv::Mat &translation_vector,
std::vector<cv::Point3d> *triangulated_points)
{
if (matches.empty())
{
std::cout << "警告: 输入匹配为空，跳过三角测量。" << std::endl;
return;
}
triangulated_points->clear();

// 投影矩阵 P1 = K[I|0], P2 = K[R|t]
// 注意: 输入的 points 已经是归一化平面坐标，因此这里的 K 实际上被抵消了。
cv::Mat proj_matrix1 = (cv::Mat_<float>(3, 4) << 1, 0, 0, 0, 0, 1, 0, 0, 0,
0, 1, 0);
cv::Mat proj_matrix2 = (cv::Mat_<float>(3, 4) << rotation_matrix.at<double>(0, 0),
rotation_matrix.at<double>(0, 1),
rotation_matrix.at<double>(0, 2),
translation_vector.at<double>(0, 0),
rotation_matrix.at<double>(1, 0),
rotation_matrix.at<double>(1, 1),
rotation_matrix.at<double>(1, 2),
translation_vector.at<double>(1, 0),
rotation_matrix.at<double>(2, 0),
rotation_matrix.at<double>(2, 1),
rotation_matrix.at<double>(2, 2),
translation_vector.at<double>(2, 0));

std::vector<cv::Point2f> normalized_points1, normalized_points2;
for (const auto &match : matches)
{
normalized_points1.push_back(
PixelToNormalizedPlane(keypoints1[match.queryIdx].pt, camera_matrix));
normalized_points2.push_back(
PixelToNormalizedPlane(keypoints2[match.trainIdx].pt, camera_matrix));
}

cv::Mat points_4d_homogeneous;
cv::triangulatePoints(proj_matrix1, proj_matrix2, normalized_points1,
normalized_points2, points_4d_homogeneous);

// 转换齐次坐标到3D坐标
for (int i = 0; i < points_4d_homogeneous.cols; ++i)
{
cv::Mat x = points_4d_homogeneous.col(i);
x /= x.at<float>(3, 0); // 归一化
// 只保留深度值为正的点
if (x.at<float>(2, 0) > 0)
{
triangulated_points->push_back(cv::Point3d(
x.at<float>(0, 0), x.at<float>(1, 0), x.at<float>(2, 0)));
}
}
}

inline cv::Scalar GetDepthColor(float depth, float upper_threshold,
float lower_threshold)
{
if (depth > upper_threshold)
depth = upper_threshold;
if (depth < lower_threshold)
depth = lower_threshold;
float range = upper_threshold – lower_threshold;
float ratio = (depth – lower_threshold) / range;
// cv::Scalar 顺序为 (B, G, R)
return cv::Scalar(255 * ratio, 0, 255 * (1 – ratio));
}

void FilterMatchesByDepthConsistency(
const std::vector<cv::Point3d> &points,
const std::vector<cv::DMatch> &input_matches,
std::vector<cv::DMatch> *filtered_matches,
double depth_threshold_factor)
{
if (points.empty())
{
std::cout << "警告: 输入点云为空，无法进行深度筛选。" << std::endl;
return;
}
filtered_matches->clear();

// 提取深度值
std::vector<double> depths;
depths.reserve(points.size());
for (const auto &p : points)
{
depths.push_back(p.z);
}

// 计算深度中位数 (Median)
std::sort(depths.begin(), depths.end());
double median_depth = depths[depths.size() / 2];

// 计算中位数绝对偏差 (MAD)
std::vector<double> abs_deviations;
abs_deviations.reserve(depths.size());
for (double d : depths)
{
abs_deviations.push_back(std::fabs(d – median_depth));
}
std::sort(abs_deviations.begin(), abs_deviations.end());
double mad = abs_deviations[abs_deviations.size() / 2];

// 设定筛选边界
double lower_bound = median_depth – depth_threshold_factor * mad;
double upper_bound = median_depth + depth_threshold_factor * mad;

// 筛选匹配点
for (size_t i = 0; i < points.size(); ++i)
{
double depth = points[i].z;
if (depth >= lower_bound && depth <= upper_bound)
{
filtered_matches->push_back(input_matches[i]);
}
}
std::cout << "深度一致性筛选: " << input_matches.size() << " -> "
<< filtered_matches->size() << " 个匹配点" << std::endl;
}

cv::Point3d FindFarthestPointPairRansac(
const std::vector<cv::Point3d> &points,
int *best_index1,
int *best_index2,
int iterations,
double inlier_threshold)
{
int n = static_cast<int>(points.size());
if (n < 2)
{
*best_index1 = *best_index2 = -1;
return cv::Point3d(0, 0, 0);
}

int best_inlier_count = 0;
double best_pair_distance = 0.0;
int local_best_idx1 = 0, local_best_idx2 = 0;

std::srand(static_cast<unsigned int>(time(nullptr)));

for (int iter = 0; iter < iterations; ++iter)
{
int i = std::rand() % n;
int j = std::rand() % n;
if (i == j)
continue;

const cv::Point3d &p1 = points[i];
const cv::Point3d &p2 = points[j];
double dist = cv::norm(p1 – p2);

int inlier_count = 0;
for (int k = 0; k < n; ++k)
{
if (k == i || k == j)
continue;
double d1 = cv::norm(points[k] – p1);
double d2 = cv::norm(points[k] – p2);
if (std::fabs(d1 – d2) < inlier_threshold * dist)
{
++inlier_count;
}
}

if (inlier_count > best_inlier_count ||
(inlier_count == best_inlier_count && dist > best_pair_distance))
{
best_inlier_count = inlier_count;
best_pair_distance = dist;
local_best_idx1 = i;
local_best_idx2 = j;
}
}

*best_index1 = local_best_idx1;
*best_index2 = local_best_idx2;
std::cout << "RANSAC找到最远点对，内点数: " << best_inlier_count << "/" << n
<< std::endl;
return cv::Point3d(points[local_best_idx1].x, points[local_best_idx1].y,
best_pair_distance);
}

bool ValidateSphereAssumption(const std::vector<cv::Point3d> &points_km)
{
if (points_km.empty())
{
return false;
}

// 计算点云质心
cv::Point3d centroid(0, 0, 0);
for (const auto &p : points_km)
{
centroid.x += p.x;
centroid.y += p.y;
centroid.z += p.z;
}
centroid.x /= points_km.size();
centroid.y /= points_km.size();
centroid.z /= points_km.size();

// 计算各点到质心的距离
double sum_radius = 0.0;
double min_radius = std::numeric_limits<double>::max();
double max_radius = 0.0;
for (const auto &p : points_km)
{
double radius = cv::norm(p – centroid);
sum_radius += radius;
if (radius < min_radius)
min_radius = radius;
if (radius > max_radius)
max_radius = radius;
}
double avg_radius = sum_radius / points_km.size();
double radius_variation = (max_radius – min_radius) / avg_radius;

std::cout << "\\n— 球面假设验证 —" << std::endl;
std::cout << "平均半径: " << avg_radius << " km" << std::endl;
std::cout << "半径范围: [" << min_radius << ", " << max_radius << "] km"
<< std::endl;
std::cout << "半径变化率: " << radius_variation * 100 << "%" << std::endl;
std::cout << "实际月球半径: " << constants::kMoonRadiusKm << " km" << std::endl;
std::cout << "平均半径误差: "
<< std::fabs(avg_radius – constants::kMoonRadiusKm) /
constants::kMoonRadiusKm * 100
<< "%" << std::endl;

const double kVariationThreshold = 0.2; // 20%
return radius_variation < kVariationThreshold;
}

int ORB_pattern[256 * 4] = {
//这里省略
};

六、运行

6.1  输出

========== 视频文件信息 ==========
文件路径：../moon_video.mp4
—————————————-
总帧数：434 帧
帧率：16.0371 FPS
分辨率：2592 x 1944 像素
宽高比：1.33333
—————————————-
编码器 (FOURCC)：828601953 (avc1)
比特率：26512 bps
格式：0
模式：0
—————————————-
亮度：0
对比度：0
饱和度：0
色调：0
—————————————-
预计播放时长：27.062244 秒
========================================

内参矩阵:
[227208.9761570828, 0, 1296;
0, 227102.8037383177, 972;
0, 0, 1]

========== 视频分半处理配置 ==========
总帧数: 434 帧
帧间隔: 2 帧
前半部分: 0 ~ 431 帧
后半部分: 2 ~ 433 帧
每对处理帧之间的时间差: 0.124711 秒
========================================

============================================================================ 处理帧对 1 / 432 ============================================================================
前帧索引: 0，后帧索引: 2

========== 第一步：特征提取 ==========
无效/总关键点: 0/251
无效/总关键点: 0/144
提取ORB特征耗时 = 0.00305805 秒。

========== 第二步：特征匹配 ==========
ORB特征匹配耗时 = 0.000169268 秒。
本帧匹配对数： 31

========== 第三步：基于对极几何估计相机运动 ==========
估计的旋转矩阵 R =
[0.99999972001521, 0.0007482648856260482, 8.316374434720681e-06;
-0.0007482649476611471, 0.9999997200219275, 7.4588472919597e-06;
-8.310790912974711e-06, -7.465068055090365e-06, 0.9999999999376017]
估计的平移向量 t =
[0.007741135626734912;
0.0004223463444822948;
-0.9999699477698186]
平均对极约束误差:4.74887e-07

========== 第四步：基于三角化测量目标点 ==========
深度一致性筛选: 31 -> 7 个匹配点
匹配点数量：31筛选后匹配点数量：7
— 归一化尺度统计 —
最小深度: 2.45215
最大深度: 12451.5
平均深度: 1568.66

========== 第五步：重建尺度 ==========
RANSAC找到最远点对，内点数: 2/10
points_norm:10idx1:6idx2:1
diameter_norm:5.11169
moon_diameter_norm:0.0103525
尺度因子1:680.065
尺度因子2:0.0502585
尺度因子3:335793
尺度因子:335793
[2599.418493430073;
141.8209099306174;
-335782.8231467805]
相机间物理距离：335793 km
推算速度：2.69257e+06 km/s
physical_distance:825151
best_distance:0,100%
重建的深度/基线比：-nan
真实的深度/基线比：7.688e+06
误差倍数：-nan

— 球面拟合验证 —
平均半径: 7.56327e+08 km
最小半径: 1.10563e+08 km
最大半径: 3.65439e+09 km
半径变化率: 468.557%
实际月球半径: 1738.14 km
半径误差: 4.35135e+07%
点云分布不均匀

========== 第六步：可视化与输出 ==========

============================================================================ 处理帧对 2 / 432 ============================================================================
前帧索引: 1，后帧索引: 3

========== 第一步：特征提取 ==========
无效/总关键点: 0/141
无效/总关键点: 0/158
提取ORB特征耗时 = 0.00304711 秒。

========== 第二步：特征匹配 ==========
ORB特征匹配耗时 = 7.5632e-05 秒。
本帧匹配对数： 21

========== 第三步：基于对极几何估计相机运动 ==========
估计的旋转矩阵 R =
[-0.9999632595362404, 0.004209017409307915, -0.00746751297992468;
-0.00424779276219489, -0.9999775393952385, 0.005184299537426044;
-0.007445524448059372, 0.005215829511444477, 0.9999586787903793]
估计的平移向量 t =
[-0.003732644367462738;
0.002597217279565931;
0.9999896608607654]
平均对极约束误差:1.94261e-07

========== 第四步：基于三角化测量目标点 ==========
深度一致性筛选: 21 -> 5 个匹配点
匹配点数量：21筛选后匹配点数量：5
— 归一化尺度统计 —
最小深度: 0.0524616
最大深度: 0.583547
平均深度: 0.361048

========== 第五步：重建尺度 ==========
RANSAC找到最远点对，内点数: 1/6
points_norm:6idx1:5idx2:1
diameter_norm:0.531091
moon_diameter_norm:0.0103525
尺度因子1:6545.55
尺度因子2:0.0502585
尺度因子3:335793
尺度因子:335793
[-1253.395530840568;
872.1271598244039;
335789.4426630427]
相机间物理距离：335793 km
推算速度：2.69257e+06 km/s
physical_distance:19354.4
best_distance:19354.4,94.965%
重建的深度/基线比：0.0576378
真实的深度/基线比：7.688e+06
误差倍数：1.33385e+08

— 球面拟合验证 —
平均半径: 39578.6 km
最小半径: 519.407 km
最大半径: 103622 km
半径变化率: 260.501%
实际月球半径: 1738.14 km
半径误差: 2177.07%
点云分布不均匀

========== 第六步：可视化与输出 ==========

（省略）

============================================================================ 处理帧对 431 / 432 ============================================================================
前帧索引: 430，后帧索引: 432

========== 第一步：特征提取 ==========
无效/总关键点: 5/782
无效/总关键点: 4/724
提取ORB特征耗时 = 0.00733437 秒。

========== 第二步：特征匹配 ==========
ORB特征匹配耗时 = 0.00147615 秒。
本帧匹配对数： 98

========== 第三步：基于对极几何估计相机运动 ==========
估计的旋转矩阵 R =
[0.9999999880755438, 0.0001541835737123092, 8.737152105603918e-06;
-0.0001541836488287737, 0.999999988076744, 8.597353473714915e-06;
-8.735826430745827e-06, -8.598700497188932e-06, 0.9999999999248738]
估计的平移向量 t =
[0.001556202400334755;
0.001404142764049561;
0.9999978033061809]
平均对极约束误差:5.33428e-07

========== 第四步：基于三角化测量目标点 ==========
深度一致性筛选: 98 -> 35 个匹配点
匹配点数量：98筛选后匹配点数量：35
— 归一化尺度统计 —
最小深度: 0.131045
最大深度: 91634
平均深度: 4651.98

========== 第五步：重建尺度 ==========
RANSAC找到最远点对，内点数: 6/43
points_norm:43idx1:34idx2:7
diameter_norm:703.262
moon_diameter_norm:0.0103525
尺度因子1:4.94308
尺度因子2:0.0502585
尺度因子3:335793
尺度因子:335793
[522.5617395178803;
471.5011910760103;
335792.1768385198]
相机间物理距离：335793 km
推算速度：2.69257e+06 km/s
physical_distance:45742.2
best_distance:384580,0.0469319%
重建的深度/基线比：1.14529
真实的深度/基线比：7.688e+06
误差倍数：6.71271e+06

— 球面拟合验证 —
平均半径: 2.32223e+09 km
最小半径: 1.19126e+09 km
最大半径: 2.92081e+10 km
半径变化率: 1206.46%
实际月球半径: 1738.14 km
半径误差: 1.33604e+08%
点云分布不均匀

========== 第六步：可视化与输出 ==========

============================================================================ 处理帧对 432 / 432 ============================================================================
前帧索引: 431，后帧索引: 433

========== 第一步：特征提取 ==========
无效/总关键点: 4/713
无效/总关键点: 3/702
提取ORB特征耗时 = 0.00742921 秒。

========== 第二步：特征匹配 ==========
ORB特征匹配耗时 = 0.00132217 秒。
本帧匹配对数： 114

========== 第三步：基于对极几何估计相机运动 ==========
估计的旋转矩阵 R =
[0.999999995693401, -9.221188040359199e-05, 1.049604853560725e-05;
9.221179067177388e-05, 0.9999999957119479, 8.549271732712782e-06;
-1.049683683502255e-05, -8.548303836464647e-06, 0.9999999999083714]
估计的平移向量 t =
[-0.004435389846736783;
0.00156595639178756;
0.9999889374875536]
平均对极约束误差:2.31351e-07

========== 第四步：基于三角化测量目标点 ==========
深度一致性筛选: 114 -> 61 个匹配点
匹配点数量：114筛选后匹配点数量：61
— 归一化尺度统计 —
最小深度: 0.321582
最大深度: 5417.31
平均深度: 898.267

========== 第五步：重建尺度 ==========
RANSAC找到最远点对，内点数: 11/91
points_norm:91idx1:9idx2:30
diameter_norm:324.903
moon_diameter_norm:0.0103525
尺度因子1:10.6994
尺度因子2:0.0502585
尺度因子3:335793
尺度因子:335793
[-1489.372483458542;
525.8370607355598;
335789.1997594438]
相机间物理距离：335793 km
推算速度：2.69257e+06 km/s
physical_distance:109723
best_distance:384580,0.0469319%
重建的深度/基线比：1.14529
真实的深度/基线比：7.688e+06
误差倍数：6.71271e+06

— 球面拟合验证 —
平均半径: 3.52298e+08 km
最小半径: 3.24829e+07 km
最大半径: 1.51748e+09 km
半径变化率: 421.517%
实际月球半径: 1738.14 km
半径误差: 2.02686e+07%
点云分布不均匀

========== 第六步：可视化与输出 ==========
后视频已读取完毕。

========== 处理摘要 ==========
处理的总帧数: 431
找到的总匹配对数: 16066
平均每帧匹配对数: 37

按任意键关闭窗口…

6.2 运行画面

可以直观看出重建不可靠

七、SLAM技术路线对比

7.1核心选择指南与趋势

7.2融合才是未来：主流方案已非纯技术路线

在实际前沿应用中，纯单目、双目或激光雷达的SLAM系统已较少见，取而代之的是多传感器融合方案，以克服单一传感器的固有缺陷：

• 视觉-惯性里程计（VIO）：单目/双目 + IMU。成为移动机器人和AR/VR的事实标准。IMU弥补了视觉的短板，提供尺度、重力参考并应对快速运动。
• 激光雷达-惯性里程计（LIO）：激光雷达 + IMU。IMU校正点云运动畸变，实现更精准、鲁棒的建图。
• 激光雷达-视觉-惯性融合（LVI-SLAM）：激光雷达 + 相机 + IMU。目前的“天花板”配置，实现最大程度的鲁棒性和精度。

没有“最好”的SLAM路线，只有“最合适”的。计算量上，双目将负担前置（立体匹配），单目将负担后置（全局优化），激光雷达则用于处理确定性的几何数据。性能上，激光雷达在精度和鲁棒性上全面领先，但代价是成本和数据处理量。