Cartographer源码阅读3D-回环计算-分支定界

计算回环-ConstraintBuilder3D

ConstraintBuilder2D基本流程一致,通过2个接口:MaybeAddConstraint和MaybeAddGlobalConstraint,然后调用统一接口:ComputeConstraint函数,计算分支定界,然后利用最高得分的解作为初始解,进行CSM位姿解算。

和2D不同的地方主要包括:

(1)分支定界的方法不2D不同:

多分辨率网格的形式不一致;

离散方式不同:先通过直方图匹配的方法计算进行角度离散,再进行距离离散;

分支操作:分支时,考虑深度,每次分支需要计算的子节点变为8个。

(2)CSM需要考虑高低分辨率的匹配。

CSM不同的地方主要体现在:

void ConstraintBuilder3D::ComputeConstraint(
    const SubmapId& submap_id, const NodeId& node_id, bool match_full_submap,
    const TrajectoryNode::Data* const constant_data,
    const transform::Rigid3d& global_node_pose,
    const transform::Rigid3d& global_submap_pose,
    const SubmapScanMatcher& submap_scan_matcher,
    std::unique_ptr<Constraint>* constraint) {
  // The 'constraint_transform' (submap i <- node j) is computed from:
  // - a 'high_resolution_point_cloud' in node j and
  // - the initial guess 'initial_pose' (submap i <- node j).
  std::unique_ptr<scan_matching::FastCorrelativeScanMatcher3D::Result>
      match_result;

  // Compute 'pose_estimate' in three stages:
  // 1. Fast estimate using the fast correlative scan matcher.
  // 2. Prune if the score is too low.
  // 3. Refine.
  if (match_full_submap) {
    kGlobalConstraintsSearchedMetric->Increment();
    // 分支定界搜索整个submap,传进去的是submap和node的global rotation及global_localization_min_score
    match_result =
        submap_scan_matcher.fast_correlative_scan_matcher->MatchFullSubmap(
            global_node_pose.rotation(), global_submap_pose.rotation(),
            *constant_data, options_.global_localization_min_score());
    if (match_result != nullptr) {
      CHECK_GT(match_result->score, options_.global_localization_min_score());
      CHECK_GE(node_id.trajectory_id, 0);
      CHECK_GE(submap_id.trajectory_id, 0);
      kGlobalConstraintsFoundMetric->Increment();
      kGlobalConstraintScoresMetric->Observe(match_result->score);
      kGlobalConstraintRotationalScoresMetric->Observe(
          match_result->rotational_score);
      kGlobalConstraintLowResolutionScoresMetric->Observe(
          match_result->low_resolution_score);
    } else {
      return;
    }
  } else {
    kConstraintsSearchedMetric->Increment();
    // 分支定界搜索部分submap,传进去的是submap和node的global pose,分数为min_score
    match_result = submap_scan_matcher.fast_correlative_scan_matcher->Match(
        global_node_pose, global_submap_pose, *constant_data,
        options_.min_score());
    if (match_result != nullptr) {
      // We've reported a successful local match.
      CHECK_GT(match_result->score, options_.min_score());
      kConstraintsFoundMetric->Increment();
      kConstraintScoresMetric->Observe(match_result->score);
      kConstraintRotationalScoresMetric->Observe(
          match_result->rotational_score);
      kConstraintLowResolutionScoresMetric->Observe(
          match_result->low_resolution_score);
    } else {
      return;
    }
  }
  {
    common::MutexLocker locker(&mutex_);
    score_histogram_.Add(match_result->score);
    rotational_score_histogram_.Add(match_result->rotational_score);
    low_resolution_score_histogram_.Add(match_result->low_resolution_score);
  }

  // Use the CSM estimate as both the initial and previous pose. This has the
  // effect that, in the absence of better information, we prefer the original
  // CSM estimate.
  // CSM计算约束,分为高低分辨率的网格和前端计算CSM基本一致
  ceres::Solver::Summary unused_summary;
  transform::Rigid3d constraint_transform;
  ceres_scan_matcher_.Match(match_result->pose_estimate.translation(),
                            match_result->pose_estimate,
                            {{&constant_data->high_resolution_point_cloud,
                              submap_scan_matcher.high_resolution_hybrid_grid},
                             {&constant_data->low_resolution_point_cloud,
                              submap_scan_matcher.low_resolution_hybrid_grid}},
                            &constraint_transform, &unused_summary);

  constraint->reset(new Constraint{
      submap_id,
      node_id,
      {constraint_transform, options_.loop_closure_translation_weight(),
       options_.loop_closure_rotation_weight()},
      Constraint::INTER_SUBMAP});

  if (options_.log_matches()) {
    std::ostringstream info;
    info << "Node " << node_id << " with "
         << constant_data->high_resolution_point_cloud.size()
         << " points on submap " << submap_id << std::fixed;
    if (match_full_submap) {
      info << " matches";
    } else {
      // Compute the difference between (submap i <- node j) according to loop
      // closure ('constraint_transform') and according to global SLAM state.
      const transform::Rigid3d difference = global_node_pose.inverse() *
                                            global_submap_pose *
                                            constraint_transform;
      info << " differs by translation " << std::setprecision(2)
           << difference.translation().norm() << " rotation "
           << std::setprecision(3) << transform::GetAngle(difference);
    }
    info << " with score " << std::setprecision(1) << 100. * match_result->score
         << "%.";
    LOG(INFO) << info.str();
  }
}

分支定界-FastCorrelativeScanMatcher3D:

和2D一样,有2个入口:MatchFullSubmap和Match,MatchFullSubmap利用默认参数,设置一个较大的窗口,而Match则是构造一个局部的小窗口,且都在最粗糙的网格进行暴力匹配,计算得分。统一调用MatchWithSearchParameters函数,先进行点云离散,生成候选解(通过直方图匹配进行角度离散,再进行xyz方向距离离散),然后进行分支定界。

class PrecomputationGridStack3D {
 public:
  PrecomputationGridStack3D(
      const HybridGrid& hybrid_grid,
      const proto::FastCorrelativeScanMatcherOptions3D& options);

  const PrecomputationGrid3D& Get(int depth) const {
    return precomputation_grids_.at(depth);
  }

  int max_depth() const { return precomputation_grids_.size() - 1; }

 private:
  // 存放由精到糙的网格
  std::vector<PrecomputationGrid3D> precomputation_grids_;
};

struct DiscreteScan3D;
struct Candidate3D;

struct DiscreteScan3D {
  transform::Rigid3f pose;
  // Contains a vector of discretized scans for each 'depth'.
  std::vector<std::vector<Eigen::Array3i>> cell_indices_per_depth;
  float rotational_score;
};

struct Candidate3D {
  Candidate3D(const int scan_index, const Eigen::Array3i& offset)
      : scan_index(scan_index), offset(offset) {}

  static Candidate3D Unsuccessful() {
    return Candidate3D(0, Eigen::Array3i::Zero());
  }

  // Index into the discrete scans vectors.
  int scan_index;

  // Linear offset from the initial pose in cell indices. For lower resolution
  // candidates this is the lowest offset of the 2^depth x 2^depth x 2^depth
  // block of possibilities.
  Eigen::Array3i offset;

  // Score, higher is better.
  float score = -std::numeric_limits<float>::infinity();

  // Score of the low resolution matcher.
  float low_resolution_score = 0.f;

  bool operator<(const Candidate3D& other) const { return score < other.score; }
  bool operator>(const Candidate3D& other) const { return score > other.score; }
};

// Used to compute scores between 0 and 1 how well the given pose matches.
using MatchingFunction = std::function<float(const transform::Rigid3f&)>;

class FastCorrelativeScanMatcher3D {
 public:
  struct Result {
    float score;
    transform::Rigid3d pose_estimate;
    float rotational_score;
    float low_resolution_score;
  };

  FastCorrelativeScanMatcher3D(
      const HybridGrid& hybrid_grid,
      const HybridGrid* low_resolution_hybrid_grid,
      const std::vector<TrajectoryNode>& nodes,
      const proto::FastCorrelativeScanMatcherOptions3D& options);
  ~FastCorrelativeScanMatcher3D();

  FastCorrelativeScanMatcher3D(const FastCorrelativeScanMatcher3D&) = delete;
  FastCorrelativeScanMatcher3D& operator=(const FastCorrelativeScanMatcher3D&) =
      delete;

  // Aligns the node with the given 'constant_data' within the 'hybrid_grid'
  // given 'global_node_pose' and 'global_submap_pose'. 'Result' is only
  // returned if a score above 'min_score' (excluding equality) is possible.
  std::unique_ptr<Result> Match(const transform::Rigid3d& global_node_pose,
                                const transform::Rigid3d& global_submap_pose,
                                const TrajectoryNode::Data& constant_data,
                                float min_score) const;

  // Aligns the node with the given 'constant_data' within the 'hybrid_grid'
  // given rotations which are expected to be approximately gravity aligned.
  // 'Result' is only returned if a score above 'min_score' (excluding equality)
  // is possible.
  std::unique_ptr<Result> MatchFullSubmap(
      const Eigen::Quaterniond& global_node_rotation,
      const Eigen::Quaterniond& global_submap_rotation,
      const TrajectoryNode::Data& constant_data, float min_score) const;

 private:
  struct SearchParameters {
    const int linear_xy_window_size;     // voxels
    const int linear_z_window_size;      // voxels
    const double angular_search_window;  // radians
    const MatchingFunction* const low_resolution_matcher;
  };

  std::unique_ptr<Result> MatchWithSearchParameters(
      const SearchParameters& search_parameters,
      const transform::Rigid3f& global_node_pose,
      const transform::Rigid3f& global_submap_pose,
      const sensor::PointCloud& point_cloud,
      const Eigen::VectorXf& rotational_scan_matcher_histogram,
      const Eigen::Quaterniond& gravity_alignment, float min_score) const;
  DiscreteScan3D DiscretizeScan(const SearchParameters& search_parameters,
                                const sensor::PointCloud& point_cloud,
                                const transform::Rigid3f& pose,
                                float rotational_score) const;
  std::vector<DiscreteScan3D> GenerateDiscreteScans(
      const SearchParameters& search_parameters,
      const sensor::PointCloud& point_cloud,
      const Eigen::VectorXf& rotational_scan_matcher_histogram,
      const Eigen::Quaterniond& gravity_alignment,
      const transform::Rigid3f& global_node_pose,
      const transform::Rigid3f& global_submap_pose) const;
  std::vector<Candidate3D> GenerateLowestResolutionCandidates(
      const SearchParameters& search_parameters, int num_discrete_scans) const;
  void ScoreCandidates(int depth,
                       const std::vector<DiscreteScan3D>& discrete_scans,
                       std::vector<Candidate3D>* const candidates) const;
  std::vector<Candidate3D> ComputeLowestResolutionCandidates(
      const SearchParameters& search_parameters,
      const std::vector<DiscreteScan3D>& discrete_scans) const;
  Candidate3D BranchAndBound(const SearchParameters& search_parameters,
                             const std::vector<DiscreteScan3D>& discrete_scans,
                             const std::vector<Candidate3D>& candidates,
                             int candidate_depth, float min_score) const;
  transform::Rigid3f GetPoseFromCandidate(
      const std::vector<DiscreteScan3D>& discrete_scans,
      const Candidate3D& candidate) const;

  const proto::FastCorrelativeScanMatcherOptions3D options_;
  const float resolution_;
  const int width_in_voxels_;
  std::unique_ptr<PrecomputationGridStack3D> precomputation_grid_stack_;
  const HybridGrid* const low_resolution_hybrid_grid_;
  RotationalScanMatcher rotational_scan_matcher_;
};

// 构造函数,构造多分辨率网格管理器:precomputation_grid_stack_,构造直方图匹配器rotational_scan_matcher
FastCorrelativeScanMatcher3D::FastCorrelativeScanMatcher3D(
    const HybridGrid& hybrid_grid,
    const HybridGrid* const low_resolution_hybrid_grid,
    const std::vector<TrajectoryNode>& nodes,
    const proto::FastCorrelativeScanMatcherOptions3D& options)
    : options_(options),
      resolution_(hybrid_grid.resolution()),
      width_in_voxels_(hybrid_grid.grid_size()),
      precomputation_grid_stack_(
          common::make_unique<PrecomputationGridStack3D>(hybrid_grid, options)),
      low_resolution_hybrid_grid_(low_resolution_hybrid_grid),
      rotational_scan_matcher_(HistogramsAtAnglesFromNodes(nodes)) {}

// 入口函数
std::unique_ptr<FastCorrelativeScanMatcher3D::Result>
FastCorrelativeScanMatcher3D::Match(
    const transform::Rigid3d& global_node_pose,
    const transform::Rigid3d& global_submap_pose,
    const TrajectoryNode::Data& constant_data, const float min_score) const {
  // 构造函数,使用submap的低分辨率网格和node的低分辨率点云,传入pose,暴力计算该pose下,点云和网格的匹配得分
  const auto low_resolution_matcher = scan_matching::CreateLowResolutionMatcher(
      low_resolution_hybrid_grid_, &constant_data.low_resolution_point_cloud);
  // 构造搜索框
  const SearchParameters search_parameters{
      common::RoundToInt(options_.linear_xy_search_window() / resolution_),
      common::RoundToInt(options_.linear_z_search_window() / resolution_),
      options_.angular_search_window(), &low_resolution_matcher};
  // 分支定界搜索
  return MatchWithSearchParameters(
      search_parameters, global_node_pose.cast<float>(),
      global_submap_pose.cast<float>(),
      constant_data.high_resolution_point_cloud,
      constant_data.rotational_scan_matcher_histogram,
      constant_data.gravity_alignment, min_score);
}
// 入口函数
std::unique_ptr<FastCorrelativeScanMatcher3D::Result>
FastCorrelativeScanMatcher3D::MatchFullSubmap(
    const Eigen::Quaterniond& global_node_rotation,
    const Eigen::Quaterniond& global_submap_rotation,
    const TrajectoryNode::Data& constant_data, const float min_score) const {
  // 计算最远点
  float max_point_distance = 0.f;
  for (const Eigen::Vector3f& point :
       constant_data.high_resolution_point_cloud) {
    max_point_distance = std::max(max_point_distance, point.norm());
  }
  // 构造最大线搜索框
  const int linear_window_size =
      (width_in_voxels_ + 1) / 2 +
      common::RoundToInt(max_point_distance / resolution_ + 0.5f);
  // 生成构造函数,使用submap的低分辨率网格和node的低分辨率点云,传入pose,暴力计算该pose下,点云和网格的匹配得分
  const auto low_resolution_matcher = scan_matching::CreateLowResolutionMatcher(
      low_resolution_hybrid_grid_, &constant_data.low_resolution_point_cloud);
  // 构造搜索框
  const SearchParameters search_parameters{
      linear_window_size, linear_window_size, M_PI, &low_resolution_matcher};
  return MatchWithSearchParameters(
      search_parameters,
      transform::Rigid3f::Rotation(global_node_rotation.cast<float>()),
      transform::Rigid3f::Rotation(global_submap_rotation.cast<float>()),
      constant_data.high_resolution_point_cloud,
      constant_data.rotational_scan_matcher_histogram,
      constant_data.gravity_alignment, min_score);
}
// 返回函数
std::function<float(const transform::Rigid3f&)> CreateLowResolutionMatcher(
    const HybridGrid* low_resolution_grid, const sensor::PointCloud* points) {
  return [=](const transform::Rigid3f& pose) {
    // 计算平均分
    float score = 0.f;
    for (const Eigen::Vector3f& point :
         // 计算位姿转换后的点云
         sensor::TransformPointCloud(*points, pose)) {
      // TODO(zhengj, whess): Interpolate the Grid to get better score.
      // 如上注释所述,应该用插值函数来计算概率,不应该直接取周围网格的概率
      score += low_resolution_grid->GetProbability(
          low_resolution_grid->GetCellIndex(point));
    }
    // 暴力得分没有加上权重(和前端的暴力匹配不同)
    return score / points->size();
  };
}
// 分支定界入口
std::unique_ptr<FastCorrelativeScanMatcher3D::Result>
FastCorrelativeScanMatcher3D::MatchWithSearchParameters(
    const FastCorrelativeScanMatcher3D::SearchParameters& search_parameters,
    const transform::Rigid3f& global_node_pose,
    const transform::Rigid3f& global_submap_pose,
    const sensor::PointCloud& point_cloud,
    const Eigen::VectorXf& rotational_scan_matcher_histogram,
    const Eigen::Quaterniond& gravity_alignment, const float min_score) const {
  // 离散点云
  const std::vector<DiscreteScan3D> discrete_scans = GenerateDiscreteScans(
      search_parameters, point_cloud, rotational_scan_matcher_histogram,
      gravity_alignment, global_node_pose, global_submap_pose);
  // 生成低分辨率候选解(打分过的)
  const std::vector<Candidate3D> lowest_resolution_candidates =
      ComputeLowestResolutionCandidates(search_parameters, discrete_scans);
  // 分支定界计算最优解
  const Candidate3D best_candidate = BranchAndBound(
      search_parameters, discrete_scans, lowest_resolution_candidates,
      precomputation_grid_stack_->max_depth(), min_score);
  if (best_candidate.score > min_score) {
    return common::make_unique<Result>(Result{
        best_candidate.score,
        GetPoseFromCandidate(discrete_scans, best_candidate).cast<double>(),
        discrete_scans[best_candidate.scan_index].rotational_score,
        best_candidate.low_resolution_score});
  }
  return nullptr;
}

std::vector<DiscreteScan3D> FastCorrelativeScanMatcher3D::GenerateDiscreteScans(
    const FastCorrelativeScanMatcher3D::SearchParameters& search_parameters,
    const sensor::PointCloud& point_cloud,
    const Eigen::VectorXf& rotational_scan_matcher_histogram,
    const Eigen::Quaterniond& gravity_alignment,
    const transform::Rigid3f& global_node_pose,
    const transform::Rigid3f& global_submap_pose) const {
  std::vector<DiscreteScan3D> result;
  // We set this value to something on the order of resolution to make sure that
  // the std::acos() below is defined.
  // 计算最远点
  float max_scan_range = 3.f * resolution_; // 防止max_scan_range太小
  for (const Eigen::Vector3f& point : point_cloud) {
    const float range = point.norm();
    max_scan_range = std::max(range, max_scan_range);
  }
  // 和2D一样,计算角度步长
  const float kSafetyMargin = 1.f - 1e-2f;
  const float angular_step_size =
      kSafetyMargin * std::acos(1.f - common::Pow2(resolution_) /
                                          (2.f * common::Pow2(max_scan_range)));
  const int angular_window_size = common::RoundToInt(
      search_parameters.angular_search_window / angular_step_size);
  // 按照步长,分割yaw角(-max,max)
  std::vector<float> angles;
  for (int rz = -angular_window_size; rz <= angular_window_size; ++rz) {
    angles.push_back(rz * angular_step_size);
  }
  // 计算node和submap之间的相对位姿
  const transform::Rigid3f node_to_submap =
      global_submap_pose.inverse() * global_node_pose;
  // 直方图匹配,通过不断旋转直方图,得到各个角度范围下直方图的匹配得分
  const std::vector<float> scores = rotational_scan_matcher_.Match(
      rotational_scan_matcher_histogram,
      // 要求是重力对齐的,所以需要再去掉重力方向
      transform::GetYaw(node_to_submap.rotation() *
                        gravity_alignment.inverse().cast<float>()),
      angles);
  for (size_t i = 0; i != angles.size(); ++i) {
    // 小于得分的,不是候选解,大于得分的为候选解
    if (scores[i] < options_.min_rotational_score()) {
      continue;
    }
    
    // 角度增量
    const Eigen::Vector3f angle_axis(0.f, 0.f, angles[i]);
    // It's important to apply the 'angle_axis' rotation between the translation
    // and rotation of the 'initial_pose', so that the rotation is around the
    // origin of the range data, and yaw is in map frame.
    // 约束加上角度增量
    const transform::Rigid3f pose(
        node_to_submap.translation(),
        global_submap_pose.rotation().inverse() *
            transform::AngleAxisVectorToRotationQuaternion(angle_axis) *
            global_node_pose.rotation());
    result.push_back(
        // 点云在该pose下离散
        DiscretizeScan(search_parameters, point_cloud, pose, scores[i]));
  }
  return result;
}

DiscreteScan3D FastCorrelativeScanMatcher3D::DiscretizeScan(
    const FastCorrelativeScanMatcher3D::SearchParameters& search_parameters,
    const sensor::PointCloud& point_cloud, const transform::Rigid3f& pose,
    const float rotational_score) const {
  std::vector<std::vector<Eigen::Array3i>> cell_indices_per_depth;
  // 获取原始网格
  const PrecomputationGrid3D& original_grid =
      precomputation_grid_stack_->Get(0);
  std::vector<Eigen::Array3i> full_resolution_cell_indices;
  for (const Eigen::Vector3f& point :
       // 点云按照pose转换
       sensor::TransformPointCloud(point_cloud, pose)) {
    // 在原始网格下每个点的cell index
    full_resolution_cell_indices.push_back(original_grid.GetCellIndex(point));
  }
  const int full_resolution_depth = std::min(options_.full_resolution_depth(),
                                             options_.branch_and_bound_depth());
  CHECK_GE(full_resolution_depth, 1);
  // 按照full_resolution_depth深度,生成cell_indices_per_depth,和高分辨地图匹配的数量
  for (int i = 0; i != full_resolution_depth; ++i) {
    cell_indices_per_depth.push_back(full_resolution_cell_indices);
  }
  // 和低分辨率地图匹配的数量
  const int low_resolution_depth =
      options_.branch_and_bound_depth() - full_resolution_depth;
  CHECK_GE(low_resolution_depth, 0);
  // 计算起点
  const Eigen::Array3i search_window_start(
      -search_parameters.linear_xy_window_size,
      -search_parameters.linear_xy_window_size,
      -search_parameters.linear_z_window_size);
  // 低分辨率地图每一层的cell index都加进来
  for (int i = 0; i != low_resolution_depth; ++i) {
    const int reduction_exponent = i + 1;
    // 计算window的起点
    const Eigen::Array3i low_resolution_search_window_start(
        search_window_start[0] >> reduction_exponent,
        search_window_start[1] >> reduction_exponent,
        search_window_start[2] >> reduction_exponent);
    // 为该层添加一个cell index的vector
    cell_indices_per_depth.emplace_back();
    for (const Eigen::Array3i& cell_index : full_resolution_cell_indices) {
      const Eigen::Array3i cell_at_start = cell_index + search_window_start;
      // 计算低分辨率的终点
      const Eigen::Array3i low_resolution_cell_at_start(
          cell_at_start[0] >> reduction_exponent,
          cell_at_start[1] >> reduction_exponent,
          cell_at_start[2] >> reduction_exponent);
      cell_indices_per_depth.back().push_back(
          // 在低分辨率网格的cell index
          low_resolution_cell_at_start - low_resolution_search_window_start);
    }
  }
  // 返回该pose下,离散点云在各个分辨率地图下的cell index,及得分
  return DiscreteScan3D{pose, cell_indices_per_depth, rotational_score};
}

std::vector<Candidate3D>
FastCorrelativeScanMatcher3D::ComputeLowestResolutionCandidates(
    const FastCorrelativeScanMatcher3D::SearchParameters& search_parameters,
    const std::vector<DiscreteScan3D>& discrete_scans) const {
  // 生成低分辨率候选解
  std::vector<Candidate3D> lowest_resolution_candidates =
      GenerateLowestResolutionCandidates(search_parameters,
                                         discrete_scans.size());
  // 低分辨率候选解打分
  ScoreCandidates(precomputation_grid_stack_->max_depth(), discrete_scans,
                  &lowest_resolution_candidates);
  return lowest_resolution_candidates;
}

std::vector<Candidate3D>
FastCorrelativeScanMatcher3D::GenerateLowestResolutionCandidates(
    const FastCorrelativeScanMatcher3D::SearchParameters& search_parameters,
    const int num_discrete_scans) const {
  // 计算所有可能解的数量
  const int linear_step_size = 1 << precomputation_grid_stack_->max_depth();
  const int num_lowest_resolution_linear_xy_candidates =
      (2 * search_parameters.linear_xy_window_size + linear_step_size) /
      linear_step_size;
  const int num_lowest_resolution_linear_z_candidates =
      (2 * search_parameters.linear_z_window_size + linear_step_size) /
      linear_step_size;
  const int num_candidates =
      num_discrete_scans *
      common::Power(num_lowest_resolution_linear_xy_candidates, 2) *
      num_lowest_resolution_linear_z_candidates;
  std::vector<Candidate3D> candidates;
  candidates.reserve(num_candidates);
  // 生成所有的候选解
  for (int scan_index = 0; scan_index != num_discrete_scans; ++scan_index) {
    for (int z = -search_parameters.linear_z_window_size;
         z <= search_parameters.linear_z_window_size; z += linear_step_size) {
      for (int y = -search_parameters.linear_xy_window_size;
           y <= search_parameters.linear_xy_window_size;
           y += linear_step_size) {
        for (int x = -search_parameters.linear_xy_window_size;
             x <= search_parameters.linear_xy_window_size;
             x += linear_step_size) {
          candidates.emplace_back(scan_index, Eigen::Array3i(x, y, z));
        }
      }
    }
  }
  CHECK_EQ(candidates.size(), num_candidates);
  return candidates;
}
// 打分
void FastCorrelativeScanMatcher3D::ScoreCandidates(
    const int depth, const std::vector<DiscreteScan3D>& discrete_scans,
    std::vector<Candidate3D>* const candidates) const {
  const int reduction_exponent =
      std::max(0, depth - options_.full_resolution_depth() + 1);
  // 遍历所有的候选解
  for (Candidate3D& candidate : *candidates) {
    int sum = 0;
    // 取出候选帧
    const DiscreteScan3D& discrete_scan = discrete_scans[candidate.scan_index];
    // 根据深度计算Offset
    const Eigen::Array3i offset(candidate.offset[0] >> reduction_exponent,
                                candidate.offset[1] >> reduction_exponent,
                                candidate.offset[2] >> reduction_exponent);
    CHECK_LT(depth, discrete_scan.cell_indices_per_depth.size());
    for (const Eigen::Array3i& cell_index :
         // 取出该深度的cell index
         discrete_scan.cell_indices_per_depth[depth]) {
      // 计算在该深度的低分辨率网格上的位置
      const Eigen::Array3i proposed_cell_index = cell_index + offset;
      // 计算总分数
      sum += precomputation_grid_stack_->Get(depth).value(proposed_cell_index);
    }
    // 平均分为得分
    candidate.score = PrecomputationGrid3D::ToProbability(
        sum /
        static_cast<float>(discrete_scan.cell_indices_per_depth[depth].size()));
  }
  // 按照分数从小到大排序
  std::sort(candidates->begin(), candidates->end(),
            std::greater<Candidate3D>());
}

// 分支定界
Candidate3D FastCorrelativeScanMatcher3D::BranchAndBound(
    const FastCorrelativeScanMatcher3D::SearchParameters& search_parameters,
    const std::vector<DiscreteScan3D>& discrete_scans,
    const std::vector<Candidate3D>& candidates, const int candidate_depth,
    float min_score) const {
  if (candidate_depth == 0) {
    // 遍历分辨率最高的一层的所有可能解
    for (const Candidate3D& candidate : candidates) {
      // 分数小于阈值的直接跳过
      if (candidate.score <= min_score) {
        // Return if the candidate is bad because the following candidate will
        // not have better score.
        return Candidate3D::Unsuccessful();
      }
      // 计算点云和submap低分辨率网格的得分
      const float low_resolution_score =
          (*search_parameters.low_resolution_matcher)(
              GetPoseFromCandidate(discrete_scans, candidate));
      // 如果当前pose的点云和submap低分辨率网格的匹配得分大于阈值,则结束分支定界
      if (low_resolution_score >= options_.min_low_resolution_score()) {
        // We found the best candidate that passes the matching function.
        Candidate3D best_candidate = candidate;
        best_candidate.low_resolution_score = low_resolution_score;
        return best_candidate;
      }
    }

    // All candidates have good scores but none passes the matching function.
    return Candidate3D::Unsuccessful();
  }

  Candidate3D best_high_resolution_candidate = Candidate3D::Unsuccessful();
  best_high_resolution_candidate.score = min_score;
  for (const Candidate3D& candidate : candidates) {
    if (candidate.score <= min_score) {
      break;
    }
    // 分支操作
    std::vector<Candidate3D> higher_resolution_candidates;
    // // 步长减半,进行分支操作(因为网格采样的时候分辨率是除以2的),所有从粗糙层到下一层,其子节点应该有8个
    const int half_width = 1 << (candidate_depth - 1);
    for (int z : {0, half_width}) {
      if (candidate.offset.z() + z > search_parameters.linear_z_window_size) {
        break;
      }
      for (int y : {0, half_width}) {
        if (candidate.offset.y() + y >
            search_parameters.linear_xy_window_size) {
          break;
        }
        for (int x : {0, half_width}) {
          if (candidate.offset.x() + x >
              search_parameters.linear_xy_window_size) {
            break;
          }
          higher_resolution_candidates.emplace_back(
              candidate.scan_index, candidate.offset + Eigen::Array3i(x, y, z));
        }
      }
    }
    // depth -1 并打分
    ScoreCandidates(candidate_depth - 1, discrete_scans,
                    &higher_resolution_candidates);
    // 得到最高得分的候选解
    best_high_resolution_candidate = std::max(
        best_high_resolution_candidate,
        BranchAndBound(search_parameters, discrete_scans,
                       higher_resolution_candidates, candidate_depth - 1,
                       best_high_resolution_candidate.score));
  }
  return best_high_resolution_candidate;
}
多分辨率网格构建:

构造:

// unit8:0-255区间
class PrecomputationGrid3D : public HybridGridBase<uint8> {
 public:
  explicit PrecomputationGrid3D(const float resolution)
      : HybridGridBase<uint8>(resolution) {}

  // Maps values from [0, 255] to [kMinProbability, kMaxProbability].
  static float ToProbability(float value) {
    return kMinProbability +
           value * ((kMaxProbability - kMinProbability) / 255.f);
  }
};

PrecomputationGridStack3D::PrecomputationGridStack3D(
    const HybridGrid& hybrid_grid,
    const proto::FastCorrelativeScanMatcherOptions3D& options) {
  CHECK_GE(options.branch_and_bound_depth(), 1);
  CHECK_GE(options.full_resolution_depth(), 1);
  precomputation_grids_.reserve(options.branch_and_bound_depth());
  // 最高层的网格分辨率和submap的高精度网格分辨率一样
  // 把概率值归一化为[0,255],转化为PrecomputationGrid3D结构
  precomputation_grids_.push_back(ConvertToPrecomputationGrid(hybrid_grid));
  // 逐层降低分辨率放入管理器中
  Eigen::Array3i last_width = Eigen::Array3i::Ones();
  for (int depth = 1; depth != options.branch_and_bound_depth(); ++depth) {
    const bool half_resolution = depth >= options.full_resolution_depth();
    // 如果为true,取原始cell的概率值,否则,取分辨率一半的概率值:CellIndexAtHalfResolution(cell_index)
    const Eigen::Array3i next_width = ((1 << depth) * Eigen::Array3i::Ones());
    const int full_voxels_per_high_resolution_voxel =
        1 << std::max(0, depth - options.full_resolution_depth());
    // 计算分辨率增量
    const Eigen::Array3i shift = (next_width - last_width +
                                  (full_voxels_per_high_resolution_voxel - 1)) /
                                 full_voxels_per_high_resolution_voxel;
    precomputation_grids_.push_back(
        // 在上层网格的基础上采样
        PrecomputeGrid(precomputation_grids_.back(), half_resolution, shift));
    last_width = next_width;
  }
}

PrecomputationGrid3D ConvertToPrecomputationGrid(
    const HybridGrid& hybrid_grid) {
  PrecomputationGrid3D result(hybrid_grid.resolution());
  for (auto it = HybridGrid::Iterator(hybrid_grid); !it.Done(); it.Next()) {
    // 概率值归一化到0-255之间
    const int cell_value = common::RoundToInt(
        (ValueToProbability(it.GetValue()) - kMinProbability) *
        (255.f / (kMaxProbability - kMinProbability)));
    CHECK_GE(cell_value, 0);
    CHECK_LE(cell_value, 255);
    *result.mutable_value(it.GetCellIndex()) = cell_value;
  }
  return result;
}

PrecomputationGrid3D PrecomputeGrid(const PrecomputationGrid3D& grid,
                                    const bool half_resolution,
                                    const Eigen::Array3i& shift) {
  PrecomputationGrid3D result(grid.resolution());
  // 遍历上一层所有网格数据
  for (auto it = PrecomputationGrid3D::Iterator(grid); !it.Done(); it.Next()) {
    // xyz三个方向分别扩展shift,并且该网格的概率值,取该3个方向上扩展的所有网格的概率最大值
    for (int i = 0; i != 8; ++i) {
      // We use this value to update 8 values in the resulting grid, at
      // position (x - {0, 'shift'}, y - {0, 'shift'}, z - {0, 'shift'}).
      // If 'shift' is 2 ** (depth - 1), where depth 0 is the original grid,
      // this results in precomputation grids analogous to the 2D case.
      const Eigen::Array3i cell_index =
          it.GetCellIndex() - shift * PrecomputationGrid3D::GetOctant(i);
      // 为true,取原始网格的概率值,为false,取x/2,y/2,z/2的值
      auto* const cell_value = result.mutable_value(
          half_resolution ? CellIndexAtHalfResolution(cell_index) : cell_index);
       // 赋最大值
      *cell_value = std::max(it.GetValue(), *cell_value);
    }
  }
  return result;
}


inline int DivideByTwoRoundingTowardsNegativeInfinity(const int value) {
  return value >> 1;
}

// xyz方向分别除以2
// Computes the half resolution index corresponding to the full resolution
// 'cell_index'.
Eigen::Array3i CellIndexAtHalfResolution(const Eigen::Array3i& cell_index) {
  return Eigen::Array3i(
      DivideByTwoRoundingTowardsNegativeInfinity(cell_index[0]),
      DivideByTwoRoundingTowardsNegativeInfinity(cell_index[1]),
      DivideByTwoRoundingTowardsNegativeInfinity(cell_index[2]));
}
直方图匹配器构造及匹配:

构造:

std::vector<std::pair<Eigen::VectorXf, float>> HistogramsAtAnglesFromNodes(
    const std::vector<TrajectoryNode>& nodes) {
  std::vector<std::pair<Eigen::VectorXf, float>> histograms_at_angles;
  for (const auto& node : nodes) {
    histograms_at_angles.emplace_back(
        // node的直方图
        node.constant_data->rotational_scan_matcher_histogram,
        // node的全局yaw角,不带重力对齐的
        transform::GetYaw(
            node.global_pose *
            transform::Rigid3d::Rotation(
                node.constant_data->gravity_alignment.inverse())));
  }
  return histograms_at_angles;
}

RotationalScanMatcher::RotationalScanMatcher(
    const std::vector<std::pair<Eigen::VectorXf, float>>& histograms_at_angles)
    : histogram_(
          Eigen::VectorXf::Zero(histograms_at_angles.at(0).first.size())) {
  // submap的直方图为各个node的直方图之和
  for (const auto& histogram_at_angle : histograms_at_angles) {
    histogram_ +=
        RotateHistogram(histogram_at_angle.first, histogram_at_angle.second);
  }
}

// 根据角度插值计算新的直方图
Eigen::VectorXf RotateHistogram(const Eigen::VectorXf& histogram,
                                const float angle) {
  const float rotate_by_buckets = -angle * histogram.size() / M_PI;
  //本来是四舍五入取整,减去0.5变成向下取整;注意是负数;
  int full_buckets = common::RoundToInt(rotate_by_buckets - 0.5f);
  //总是一个小于1大于0的数
  const float fraction = rotate_by_buckets - full_buckets;
  while (full_buckets < 0) {
    full_buckets += histogram.size();
  }
  Eigen::VectorXf rotated_histogram_0 = Eigen::VectorXf::Zero(histogram.size());
  Eigen::VectorXf rotated_histogram_1 = Eigen::VectorXf::Zero(histogram.size());
  // 当前角度对应的rotate_by_buckets,它位于[full_buckets,full_buckets+1]u区间内,与下界的差值为fraction,所以,分别计算下界对应的直方图,和上界对应的直方图,然后用插值的方法,得到真实角度对应的直方图
  for (int i = 0; i != histogram.size(); ++i) {
    // 计算下界直方图
    rotated_histogram_0[i] = histogram[(i + full_buckets) % histogram.size()];
    // 计算上界直方图
    rotated_histogram_1[i] =
        histogram[(i + 1 + full_buckets) % histogram.size()];
  }
  // 插值得到真实的直方图
  return fraction * rotated_histogram_1 +
         (1.f - fraction) * rotated_histogram_0;
}

匹配:

std::vector<float> RotationalScanMatcher::Match(
    const Eigen::VectorXf& histogram, const float initial_angle,
    const std::vector<float>& angles) const {
  std::vector<float> result;
  result.reserve(angles.size());
  for (const float angle : angles) {
    // 按照增大的角度,旋转直方图
    const Eigen::VectorXf scan_histogram =
        RotateHistogram(histogram, initial_angle + angle);
    // 直方图匹配,计算每个角度和当前的submap的直方图之间的得分
    result.push_back(MatchHistograms(histogram_, scan_histogram));
  }
  return result;
}
// 计算两个直方图之间的余弦距离,得到相似性
float MatchHistograms(const Eigen::VectorXf& submap_histogram,
                      const Eigen::VectorXf& scan_histogram) {
  // We compute the dot product of normalized histograms as a measure of
  // similarity.
  // 归一化
  const float scan_histogram_norm = scan_histogram.norm();
  const float submap_histogram_norm = submap_histogram.norm(); 
  const float normalization = scan_histogram_norm * submap_histogram_norm;
  if (normalization < 1e-3f) {
    return 1.f;
  }
     // 计算点乘,得到两个直方图之间的cos(theta),即余弦距离
  return submap_histogram.dot(scan_histogram) / normalization;
}
  • 1
    点赞
  • 10
    收藏
    觉得还不错? 一键收藏
  • 2
    评论
评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值